An  Investigation  of  Momentum  Exchange  Parameterizations 


and  Atmospheric  Forcing  for  the  Coastal  Mixing  and  Optics  Program 

by 

Michiko  J.  Martin 

M.S.,  Troy  State  University  (1996) 

B.S.,  United  States  Naval  Academy  (1991) 

Submitted  in  partial  fulfillment  of  the 
requirements  for  the  degree  of 

Master  of  Science 
at  the 

Massachusetts  Institute  of  Technology 
and  the 

Woods  Hole  Oceanographic  Institution 
September  1998 

©  Michiko  J.  Martin,  1998.  All  rights  reserved. 

The  author  hereby  grants  to  MIT,  WHOI,  and  the  United  States  Navy  permission  to 
reproduce  and  to  distribute  copies  of  this  thesis  document  in  whole  or  in  part. 


MASSACHUSETTS  INSTITUTE 
_ OF  TECHNOLOGY 


JUN  21  1999 


LIBRARIES 

trig* 


Signature  of  Author. 


Joint  Program  in  Applied  Ocean  Physics  and  Engineering 
Massachusetts  Institute  of  Technology  and  Woods  Hole  Oceanographic  Institution 


Certified  by . 

U  James  B.  Edson 

Associate  Scientist,  Woods  Hole  Oceanographic  Institution 

Thesis  Supervisor 


Certified  by .  _  . 

Steven  P.  Anderson 

Associate  Scientist,  Woods  Hole  Oceanographic  Institution 

Thesis  Reader 


Approved  by .  ^ . 

Michael  S.  Triantafyllou 
Chairman,  Joint  Committee  for  Applied  Ocean  Physics  and  Engineering 
Massachusetts  Institute  of  Technology  and  Woods  Hole  Oceanographic  Institution 


Report  Documentation  Page 


Form  Approved 
OMB  No.  0704-0188 


Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 
VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 


1.  REPORT  DATE 

SEP  1998 


2.  REPORT  TYPE 


3.  DATES  COVERED 

00-00-1998  to  00-00-1998 


5a.  CONTRACT  NUMBER 


5b.  GRANT  NUMBER 


5c.  PROGRAM  ELEMENT  NUMBER 


5d.  PROJECT  NUMBER 


5e.  TASK  NUMBER 


5f.  WORK  UNIT  NUMBER 


4.  TITLE  AND  SUBTITLE 

An  Investigation  of  Momentum  Exchange  Parameterizations  and 
Atmospheric  Forcing  for  the  Coastal  Mixing  and  Optics  Program 

6.  AUTHOR(S) 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES)  8.  PERFORMING  ORGANIZATION 

Massachusetts  Institute  of  Technology, 77  Massachusetts  report  number 

Avenue, Cambridge, MA, 02139-4307 

9.  SPONSORING/MONITORING  AGENCY  NAME(S )  AND  ADDRESS(ES )  10.  SPONSOR/MONITOR' S  ACRONYM(S) 

11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 

12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

14.  ABSTRACT 

This  thesis  presents  an  investigation  of  the  influence  of  surface  waves  on  momentum  exchange.  A 
quantitative  comparison  of  direct  covariance  friction  velocity  measurements  to  bulk  aerodynamic  and 
inertial  dissipation  estimates  indicates  that  both  indirect  methods  systematically  underestimate  the 
momentum  flux  into  developing  seas.  To  account  for  wave-induced  processes  and  yield  improved  flux 
estimates,  modifications  to  the  traditional  flux  parameterizations  are  explored.  Modification  to  the  bulk 
aerodynamic  method  involves  incorporating  sea  state  dependence  into  the  roughness  length  calculation. 
For  the  inertial  dissipation  method,  a  new  parameterization  for  the  dimensionless  dissipation  rate  is 
proposed.  The  modifications  lead  to  improved  momentum  flux  estimates  for  both  methods. 

15.  SUBJECT  TERMS 


16.  SECURITY  CLASSIFICATION  OF: 

17.  LIMITATION  OF 

18.  NUMBER 

19a.  NAME  OF 

ABSTRACT 

OF  PAGES 

RESPONSIBLE  PERSON 

a.  REPORT 

unclassified 

b.  ABSTRACT 

unclassified 

c.  THIS  PAGE 

unclassified 

Same  as 
Report  (SAR) 

83 

Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


2 


An  Investigation  of  Momentum  Exchange  Parameterizations 


and  Atmospheric  Forcing  for  the  Coastal  Mixing  and  Optics  Program 

by 

Michiko  J.  Martin 

Submitted  in  partial  fulfillment  of  the  requirements  for  the  degree  of 
Master  of  Science  at  the  Massachusetts  Institute  of  Technology  and  the 
Woods  Hole  Oceanographic  Institution. 

September  1998 


Abstract 

This  thesis  presents  an  investigation  of  the  influence  of  surface  waves  on  momentum 
exchange.  A  quantitative  comparison  of  direct  covariance  friction  velocity  measurements 
to  bulk  aerodynamic  and  inertial  dissipation  estimates  indicates  that  both  indirect 
methods  systematically  underestimate  the  momentum  flux  into  developing  seas.  To 
account  for  wave-induced  processes  and  yield  improved  flux  estimates,  modifications  to 
the  traditional  flux  parameterizations  are  explored. 

Modification  to  the  bulk  aerodynamic  method  involves  incorporating  sea  state 
dependence  into  the  roughness  length  calculation.  For  the  inertial  dissipation  method,  a 
new  parameterization  for  the  dimensionless  dissipation  rate  is  proposed.  The 
modifications  lead  to  improved  momentum  flux  estimates  for  both  methods. 


Thesis  Supervisor:  Dr.  James  B.  Edson 

Title:  Associate  Scientist,  Woods  Hole  Oceanographic  Institution 


3 


Acknowledgments 


I  would  like  to  express  my  sincere  thanks  to  my  advisor  Dr.  Jim  Edson.  His 
motivation  and  enthusiasm  are  inspiring,  and  his  dedication  to  science  and  research  is 
matched  by  his  ability  and  patience  to  teach. 

I  would  also  like  to  thank  Dr.  Steve  Anderson  for  his  support  throughout  this  project 
and  his  help  as  a  reader  of  this  thesis. 

I  am  also  grateful  to  the  staff  and  faculty  of  the  Ocean  Engineering  Department  at 
MIT  and  the  Applied  Ocean  Physics  and  Engineering,  Education,  and  Physical 
Oceanography  Departments  at  WHOI,  particularly  Jean  Sucharewicz,  Julia  Westwater, 
and  Drs.  Mark  Grosenbaugh,  A1  Plueddemann,  Steve  Lentz  and  Mark  Baumgartner. 

The  love,  warmth,  and  support  of  my  family  have  sustained  and  stimulated  me.  This 
is  dedicated  to  my  parents  for  their  continual  support,  encouragement,  and  confidence. 


This  project  was  funded  by  the  Oceanographer  of  the  Navy. 


Table  of  Contents 


1 .  Introduction 

2.  Experiment 

2.1.  Deployment  and  Site  Selection 

2.2.  CMO  Instrumentation  and  Data  Processing 

2.2.1 .  Vector  Averaging  Wind  Recorder 

2.2.2.  Sonic  Anemometer 

2.2.3.  Wave  Height  Measurements 

2.3.  Marine  Boundary  Layer  Experiment 

3.  Theory 

3.1.  Basic  Equations 

3.2.  Equations  for  Turbulent  Mean  Flow 

3.3.  Equations  for  Turbulent  Kinetic  Energy 

4.  Methodology 

4. 1 .  Direct  Techniques  to  Determine  the  Momentum  Flux 

4.1.1.  Direct  Covariance  (Eddy  Correlation)  Method 

4. 1 .2.  Motion-corrected  Covariance  Method 

4.2.  Indirect  Techniques  to  Determine  the  Momentum  Flux 

4.2. 1 .  Mean  Profile  Method 

4.2.2.  Bulk  Aerodynamic  Method 

4.2.3.  Direcet  and  Inertial  Dissipation  Methods 

5 .  Analysis  and  Results 

5.1.  Data  Selection 

5.2.  Traditional  Approach 

5.3.  Modified  Inertial  Dissipation  Method 

5.4.  Modified  Bulk  Aerodynamic  Method 

5.4.1.  Scaling  Roughness  Length  with  Wave  Amplitude 

5.4.2.  Modified  Charnock’s  Relation 

6.  Conclusions 


5 


1.  Introduction 


Successful  investigations  of  air-sea  interactions  across  the  ocean  surface  depend  on 
accurate  measurements  of  the  exchange  of  momentum  and  energy  between  the  air  and 
ocean.  The  exchange  of  heat  energy  between  the  atmosphere  and  the  ocean  is  driven  by 
molecular,  turbulent,  and  radiative  processes  and  is  given  by  the  surface  heat  budget 

QsEA  ~  Q  -Qh~Qe  (1-1) 

s|t 

where  QSEA  is  the  net  heat  input  into  the  ocean;  Q  represents  the  net  downward  radiation 
(irradiance);  and,  QH  and  QE  are  the  upward  fluxes  of  sensible  and  latent  heat, 
respectively.  The  radiative  processes  at  the  surface  can  be  further  expressed  as 

Q*=SW-LW  (1.2) 

where  the  net  irradiance  Q*  results  from  the  absorption  of  solar  (shortwave)  radiation 
SW  less  the  loss  to  the  atmosphere  through  net  infrared  (longwave)  radiative  transfer 
LW  (Bradley  et  al.,  1991).  The  radiative  forcing  of  the  atmosphere  produces 
temperature  and  pressure  gradients  that  result  in  the  generation  of  winds  of  the  order  of 
10  ms'1  and  is  responsible  for  large-scale  atmospheric  circulations. 

The  exchange  of  sensible  QH  and  latent  QF  heat  occur  as  a  result  of  conduction, 
turbulence,  and  convection.  Conduction  is  important  only  within  the  molecular  sublayer, 
within  approximately  1  millimeter  of  the  ocean’s  surface,  where  fluid  motions  are 
strongly  suppressed  by  viscosity.  Immediately  above  the  molecular  boundary  layer  but 
below  the  mixed  layer,  which  is  characterized  by  deep  convection  and  thermally-driven 
circulations,  is  the  near-surface  layer  for  which  microscale  turbulence  due  to  shear  is  the 
dominant  type  of  fluid  motion.  Sensible  heat  and  water  vapor  are  transferred  via 
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microscale  eddies,  which  are  induced  mechanically  by  flow  over  irregularities  in  the 
surface  (Wallace  and  Hobbs,  1977),  and  larger  eddies  which  feed  off  the  mean  wind  shear 
(velocity  gradient).  Above  this  shear  layer,  turbulent  eddies  are  generated  by  a 
combination  of  shear  and  buoyancy.  In  turn,  these  turbulent  eddies  modify  the 
instantaneous  shear  and  give  winds  their  “gusty”  nature  by  transferring  momentum  by 
bodily  movement  of  air  parcels,  i.e.,  by  fast-moving  parcels  moving  downward  and  slow- 
moving  parcels  moving  upward  (Gill,  1 982).  Since  momentum  is  mass  times  velocity,  its 


transfer  per  unit  area  per  unit  time  can  be  expressed  as  a  flux  F  which  has  the  SI  units 


,  which  is  identical  to  — - ,  the  units  for  stress.  This  transfer  of  momentum,  or 


shear  stress,  from  the  winds  to  the  ocean  drives  the  waves  and  wind-driven  ocean  currents 
and  provides  kinetic  energy  for  mixing  in  the  ocean  mixed  layer. 

Together  with  radiative  transfer,  the  exchange  of  sensible  and  latent  heat  are 
important  mechanisms  for  redistributing  energy  between  the  atmosphere  and  ocean.  An 
understanding  of  the  dynamics  and  magnitude  of  this  energy  transfer  is  required  for  the 
initialization  and  excursion  of  coupled  air-ocean  models.  Even  the  most  crude  global 
circulation  models  must  place  a  lower  boundary  condition  on  the  mean  atmospheric 
variables  and  fluxes.  With  ever  growing  computational  power  and  sophistication  of 
numerical  models,  the  importance  of  incorporating  a  realistic  atmospheric  boundary  layer 
increases  tremendously.  Current  research  ( e.g .,  Webster  and  Lukas,  1992;  Miller  et  al., 
1992)  continues  to  demonstrate  that  more  accurate  parameterization  of  these  fluxes 
directly  leads  to  improved  model  simulations.  Whereas  an  additional  input  of  only 
1 0  Wm'2  to  the  ocean  and  troposphere  has  significant  climatological  effects,  the 
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uncertainty  in  current  estimates  of  the  surface  heat  budget  exceeds  this  value,  with 
discrepancies  as  large  as  80  Wm'2  in  the  warm  pool  region  of  the  western  Pacific  Ocean 
(Bradley  et  al.,  1991).  Decreasing  this  uncertainty  and  providing  improved  estimates  of 
the  heat  budget  will  result  in  increased  performance  in  numerical  models  and  improve  our 
ability  to  diagnose  and  predict  climate  and  climate  variability. 

The  research  presented  in  this  thesis  is  an  investigation  of  the  influence  of  wave- 
induced  processes  on  momentum  exchange.  The  thesis  explores  modifications  to 
traditional  flux  parameterizations  to  account  for  sea  state.  This  is  accomplished  by 
investigating  two  separate  approaches  to  incorporate  sea  state  dependence  into  the 
roughness  length  determination  in  the  bulk  aerodynamic  method  and  a  new  approach  to 
incorporate  wave-induced  effects  in  the  dimensionless  dissipation  rate  in  the  inertial 
dissipation  method.  Both  investigations  rely  on  comparison  of  the  parameterized  fluxes 
with  direct  covariance  fluxes. 

Details  of  the  experiment,  instrumentation,  and  data  processing  are  described  in 
Chapter  2.  Chapters  3  and  4  provide  the  theoretical  framework  upon  which  this 
investigation  is  based.  A  detailed  derivation  of  the  turbulent  kinetic  energy  budget  is 
presented  in  Chapter  3,  and  a  summary  of  the  four  standard  techniques  for  estimating  the 
momentum  flux,  including  an  introduction  to  the  well-accepted  Monin-Obukhov 
similarity  theory,  is  provided  in  Chapter  4.  Alternate  parameterizations  for  both  the  bulk 
aerodynamic  and  inertial  dissipation  methods  are  presented  in  Chapter  5,  and  their 
performance  is  evaluated.  The  final  chapter  summarizes  the  findings  of  this  thesis  and 
suggests  areas  of  future  research. 
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2.  Experiment 


2.1.  Deployment  and  Site  Selection 

The  data  used  in  this  investigation  were  collected  during  the  Coastal  Mixing  and 
Optics  (CMO)  experiment,  an  Accelerated  Research  Initiative  supported  by  the  Office  of 
Naval  Research  (Williams,  1997;  Galbraith  et  al.,  1997).  Four  moorings  were  deployed 
from  August  1996  to  June  1997  in  the  Middle  Atlantic  Bight,  90  km  south  of  Martha’s 
Vineyard  at  approximately  40.5°N,  70.5°E,  as  depicted  in  Figure  2.1 .  The  four  moorings 
included  a  central  mooring  on  the  70  m  isobath,  an  inshore  mooring  located 
approximately  10  km  inshore  along  the  60  m  isobath,  an  offshore  mooring  located  10  km 
offshore  along  the  80  m  isobath,  and  an  alongshore  mooring  located  25  km  east  of  the 
central  buoy  along  the  70  m  isobath.  The  region  has  fine  sediment  and  a  smooth  bottom 
that  slopes  gently  towards  the  shelf  break  40  km  to  the  south.  The  large-scale  surface 
heating  and  buoyancy  advection  at  this  site  create  a  strong  thermocline/pycnocline  during 
the  late-spring  and  summer  months.  The  extensive  deployment  was  chosen  so  that  the 
observation  period  included  the  destruction  of  thermal  stratification  of  the  water  column 
in  the  fall  and  the  subsequent  redevelopment  of  stratification  in  the  spring.  In  addition  to 
recording  the  seasonal  evolution  of  stratification,  the  passage  of  numerous  frontal  systems 
as  well  as  the  transits  of  Hurricanes  Edouard  and  Hortense  were  observed.  Due  to  the 
variety  of  weather  phenomena  experienced  at  the  site  during  the  eleven-month 
deployment,  this  data  set  represents  a  wide  range  of  conditions  and  is  unequaled  in 
previous  investigations  of  mixing  processes  in  coastal/shelf  waters. 
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2.2.  Coastal  Mixing  and  Optics  (CMO)  Instrumentation  and  Data  Processing 
While  all  moorings  were  outfitted  with  environmental  sensors,  the  data  used  in  this 
study  was  obtained  from  the  instrumentation  on  board  the  CMO  central  discus.  The 
central  three-meter  buoy,  shown  in  Figure  2.2,  supported  sensors  to  measure  a  variety  of 
oceanic  and  atmospheric  properties.  The  buoy  had  a  draft  of  about  1 .5  m  and  extended 
approximately  3.5  m  above  the  water  line.  The  buoy  was  outfitted  with  a  large  steering 
vane  mounted  to  one  of  the  mast  legs  that  served  to  minimize  flow  distortion  by  generally 
providing  the  best  exposure  of  the  meteorological  instruments. 

2.2. 1 .  Vector  Averaging  Wind  Recorder 

The  central  buoy  carried  two,  redundant  Vector  Averaging  Wind  Recorder  (VAWR) 
meteorological  packets  (Weller  et  al.,  1990).  The  VAWR  sensors  measured  wind  speed 
and  direction,  air  and  sea  surface  temperatures,  barometric  pressure,  relative  humidity, 
and  incoming  short-wave  and  long- wave  radiation  at  a  15 -minute  sampling  rate 
(Baumgartner  and  Anderson,  1997).  The  sensor  specifications  for  the  VAWR 
meteorological  package  are  listed  in  Table  1. 


Table  1.  VAWR  sensor  specifications.  Wind  speed  and  direction  are  vector  averaged  over  the  sampling  interval. 
Heights  are  reported  as  meters  above  the  mean  water  line.  [From  Baumgartner  and  Anderson,  1997  (Table 
5).l 


Parameter 

Sampling  Interval 

Range 

Accuracy 

Height 

Wind  speed 

15  minute  average 

0.2-50.0  ms'1 

±2%  above  0.7  ms'1 

Wind  direction 

15  minute  average 

0-360° 

±5.6° 

Air  temperature 

3.75  minute  average 

-10-35°C 

±0.2  °C  when  wind  >  5  ms'1 

2.63 

Sea  temperature 

3.75  minute  average 

-5-30  °C 

±0.005  °C 

-1.46 

Barometric  pressure 

2.5  second  sample 

0-1034  mbar 

±0.2  mbar  when  wind  <  20  ms'1 

2.72 

Relative  humidity 

3.5  second  sample 

0-100% 

±2% 

2.68 

Short-wave  radiation 

15  minute  average 

0-140  W  m'2 

±3%  of  value 

3.39 

Long-wave  radiation 

1 5  minute  average 

0-700  W  m'2 

±10% 

10 


Figure  2.1.  Coastal  Mixing  and  Optics  Site.  |From 
Stenner,  1 996.  | 


Figure  2.2.  Coastal  Mixing  and  Optics  central  3  m 
discus  buoy.  |From  Baumgartner  and  Anderson, 
1997  (Figure  4).] 
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2.2.2.  Sonic  Anemometer 


The  central  buoy  also  held  a  sonic  anemometer  which  measured  the  longitudinal, 
lateral  and  vertical  wind  velocities  as  well  as  the  speed  of  sound.  The  sonic  was  mounted 
between  the  VAWR  wind  sensors  at  a  height  3.3  m  above  the  sea  surface  and  was 
configured  to  burst  sample  twice  per  hour  for  a  duration  of  15  minutes. 

Manufactured  by  Gill  Instruments  Ltd.,  the  sonic  consists  of  a  sensing  head  and  three 
ultrasonic  transducer  pairs  supported  by  three  struts.  The  struts  are  placed  asymmetrically 
to  provide  ±  120°  of  open  exposure.  The  sonic  anemometer  is  specified  to  operate  with 
1.5%  accuracy  for  wind  speeds  between  0-60  ms'1  (Yelland  et  al.,  1994).  Each  pair  of 
transducers  alternately  transmits  and  receives  pulses  of  high  frequency  sound  waves.  All 
three  axes  complete  these  firings  every  6  ms,  which  represents  a  sampling  rate  of  168  Hz. 
Using  overlap  section  averaging,  every  eight  complete  transit  time  evaluations  are 
averaged  to  produce  the  21  Hz  output.  The  transit  time  data  were  transformed  into  an 
orthogonal  coordinate  system  by  which  the  speed  of  sound  and  wind  velocity  were 
calculated  by  determining  the  times  of  flight,  t;  and  />,  between  any  pair  of  transducers. 
Specifying  V  as  the  air  flow  parallel  to  the  line  of  the  transducer  pair,  the  flight  times  are 
defined  by 


*.= 


/ 


c  +  V 


and 


^2  “ 


l 


c-V 


(2.1) 


where  l  is  the  distance  between  the  transducers  (0.15  m)  and  c  is  the  speed  of  sound  in  air. 
The  wind  velocity  and  the  speed  of  sound  can  be  determined  explicitly: 


yj- 

2 


V^l  ^2  J 


(2.2) 
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(2.3) 


/ 

c  =  — 
2 


<1  1 
- 1 - 

V*1  ) 


Once  the  speed  of  sound  c  is  determined,  the  virtual  temperature  Tv  expressed  in 
Kelvin  can  be  calculated  by 


t=c2W 


403 


(2.4) 


where  both  c  and  wN  are  given  in  ms'1,  and  the  term 

WN2  =V2  +|(m  +  w)2  (2.5) 

corrects  for  the  curvature  in  the  path  that  the  sound  waves  travel  due  to  air  flow 
perpendicular  to  the  line  of  the  transducers  (Fairall  et  al,  1997). 

2.2.3 .  Wave  Height  Measurements 


The  surface  wave  data  were  measured  by  a  SeaTex  Waverider  buoy  (Barstow  et  al ., 
1991).  The  Waverider  buoy  was  configured  to  burst  sample  once  per  hour  at  a  1  Hz 
sampling  rate  for  a  duration  of  approximately  17  minutes  and  measured  motion,  spectral 
data,  and  scalar  wavefield  parameters,  including  significant  wave  height  and  primary 
wave  period  T  (Galbraith  et  al.,  1997).  The  wave  dispersion  equation  relates  the  radian 

271  271 

frequency,  o  =  — ,  and  wave  number,  k  =  — ,  to  the  water  depth,  D: 

T  'k 

(D2  =  gk  tanh(£D)  (2.6) 

where  X  is  the  wavelength  and  g  is  the  acceleration  due  to  gravity  (Earle  and  Bishop, 
1984).  Since  a  wave  travels  a  distance  of  one  wavelength  in  one  wave  period,  the  phase 
speed  can  be  defined  as 
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^tanh(LD) 


(2.7) 


Thus,  with  knowledge  of  only  the  water  depth  and  the  wave  period,  the  radian  frequency, 
wave  number,  and  phase  speed  can  be  determined  in  an  iterative  fashion. 

2.3.  Marine  Boundary  Layer  Experiment 

This  investigation  also  relies  on  data  analyzed  as  part  of  the  ONR's  Marine  Boundary 
Layers  (MBL)  Accelerated  Research  Initiative  (Edson  et  al.,  1997).  This  program  was 
designed  to  examine  the  3-D  structure  of  the  marine  boundary  layers.  The  open  ocean 
component  of  the  MBL  program  was  conducted  off  Monterey,  California,  in  April  and 
May  of  1995  and  involved  the  R/P  FLIP,  the  NOAA/Oak  Ridge  LongEZ  aircraft,  and 
three  research  vessels.  The  emphasis  of  this  program  was  on  improving  the 
understanding  of  flux  profile  relationships  over  the  open  ocean.  During  the  field 
experiment,  investigators  deployed  a  18  m  mast  at  the  end  of  FLIP'S  port  booms,  setting 
up  a  vertical  array  of  12  levels  of  cup/vane  anemometers  to  provide  mean  profiles.  Four 
levels  of  sonic  anemometer/thermometers  and  three  levels  of  static  pressure  sensors  were 
also  deployed  from  this  mast  in  order  to  examine  the  flux  and  dissipation  rate  profiles  and 
provide  the  necessary  scaling  parameters  to  compute  the  dimensionless  profile  functions. 
A  single  wave- wire  was  deployed  directly  beneath  the  mast  to  measure  the  instantaneous 
wave  height,  which  was  nominally  3  m  beneath  the  lowest  cup  anemometer. 
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3.  Theory 


3.1.  Basic  Equations 

Using  the  results  obtained  by  Shaw  (1990),  the  Navier-Stokes  equations  in  a  rotating 
Eulerian  reference  frame  are 


?ui  +^TP“'“*  =—^r(Phk  -<)-Pg*»  - 2pQs(^r| jUk 


dt 


dxv 


(3.1) 


where 


a*  =  Ti 


du  duk  2du,  ]  du,  _ 

T“  +  T“_7— +£  — S*  -  aik  +  pdlk 

dxk  dx,  3  ox  i  )  ox. 


(3.2) 


In  this  expression,  r\  is  the  dynamic  viscosity,  Q  is  the  angular  velocity  of  the  earth,  q . 

is  a  unit  vector  along  the  axis  of  rotation,  defines  the  stress  tensor,  and  £  is  a  second 

molecular  viscosity  which  applies  to  compressible  flow.  In  general,  the  atmosphere 
boundary  layer  may  be  treated  as  incompressible  flow  and  the  Boussinesq  approximations 
may  be  applied  (Panofsky  and  Dutton,  1 984).  This  reduces  the  continuity  equation 


dp  d 

—  +  • — -p«.  =0 
dt  dx. 


(3.3) 


to 


du , 


dX: 


!-  =  0 


(3.4) 


Applied  to  equation  (3.1),  these  approximations  yield 
du ,  du ,  10/?  d2u 

+ Uk  ~dx~ + 2Q£‘JklC]JUk  =  ~~~£r~ g5«  + v 


dt 


P  dx, 


dxkdxk 


(3.5) 


where  v  =  r\/p  is  the  kinematic  viscosity  (Shaw,  1990). 
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The  first  law  of  thermodynamics  applied  to  a  fluid  provides  the  general  equation  of 


heat  transfer 


P  T 


fds  ds  A 

—  +  Uf  — 

\dt 


,  du,  , 
-  ait  —— + 


dT  8Rl 


dxk  dxk  dxk  dxk 


(3.6) 


Rk  is  the  radiative  flux;  kc  is  the  molecular  conductivity;  and  s  is  entropy,  for  which 
changes  can  be  expressed  in  terms  of  the  standard  definitions  of  specific  heat  at  constant 
pressure  cp  and  the  equation  of  state 


Tds  =  c  dT  -—dp 


(3-7) 


A  common  expression  for  the  equation  of  state  is 


P  =  PRJV 


(3.8) 


where  Rd  is  the  gas  constant  for  dry  air  and  Tv  is  virtual  temperature,  a  function  of  specific 
humidity  q. 

Substituting  from  the  equation  of  state,  defining  the  potential  temperature  9  as 
0  s=  r(l000/ p)'^Cp  where  T  is  in  Kelvin  and  the  pressure  is  given  in  millibars,  integrating 
equation  (3.7),  and  neglecting  the  term  representing  mechanical-to-thermal  energy 
conversion  transforms  equation  (3.6)  to 


dQ  60  50 

■  + 


0  kr  d2T 


1  0  dRt 


dt  dt  dxk  T  pc  dxkdxk  pc  T  dxk 


(3.9) 


The  equation  for  the  conservation  of  specific  humidity  q  may  be  similarly  written  as 

(3.10) 


dq  dq  dq  _  dx q 

—  =  —  +  u, — -  =  Da - - 

dt  dt  dx,  dx  dx 


J  J 


where  Dq  is  the  molecular  diffusivity  of  water  vapor. 
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3.2.  Equations  for  Turbulent  Mean  Flow 


Equations  (3.3),  (3.5),  (3.8),  (3.9),  and  (3.10)  form  a  closed  set  of  seven  equations 
governing  seven  variables.  However,  the  equations  contain  nonlinear  terms  which  render 
them  impractical  to  describe  turbulent  flow.  Instead,  realizing  turbulent  flow  as  a 
stochastic  process  and  applying  Reynolds  decomposition,  which  expresses  the  variables 
as  a  sum  of  a  mean  and  fluctuating  part  [i.  e. ,  <j)(x,  y,  z,  t )  =  <j>(x,  y,  z,  t )  +  <)>'(x,  y,  z,t)], 
transforms  equation  (3.5)  to 


dut  du'  — du,  ,  du'  — du'  ,  du'  1  dp  \  dp' 

— L  +  — l  +  m, — L  +  «' — L  +  «, — —  +  u'i  ■  '  -  y  y 


dt  dt  J  dx , 


J  dxj  "J  dXj 


J  dX:  p  cbCj  p  dx. 


p'  1  dp  s  d(w(+w')  —  , 

+ - - - gS/3+v  ,  ,  -2QZyk\\jUk -2Qzijkx\jUk 


(3.11) 


P  P  fai 


dxjdxj 


Ensemble  averaging  (3.1 1)  yields  the  balance  for  mean  momentum  in  a  turbulent 
atmosphere: 


dut  —dut  —  l  dp  d  du{  duty 

— L  +  w, — L  +  2Q£iikr\iuk  =-— — +  v - - —L 

dt  1  dx  -  1  J  p  dxi  dXj  dxj  dx . 


(3.12) 


The  last  term  on  the  left-hand  side  of  (3.12)  represents  acceleration  and  the  Coriolis 
effect.  The  first  term  on  the  right-hand  side  is  the  pressure  gradient  force,  and  the  second 
term  is  the  frictional  force  due  to  molecular  viscosity.  The  last  term  represents  the  drag 
resulting  from  momentum  transfer  by  turbulence  and  is  the  kinematic  form  of  the 
divergence  of  the  Reynolds  stress  tensor, 


TV  =  ~PUiUj 


(3.13) 


For  boundary  layer  studies  the  most  important  terms  of  the  Reynolds  stress  tensor 
represent  the  vertical  transport  of  horizontal  momentum,  defined  as 
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(3.14) 


T  =  — pjl  u'w'  +  jv'w' 

Determining  the  means  to  best  parameterize  this  flux  is  the  main  focus  of  this  thesis  and 
will  be  discussed  in  detail  in  the  following  Chapter. 

The  same  Reynolds  decomposition  and  ensemble  averaging  procedure  that  were 
applied  to  equation  (3.5)  can  be  applied  to  the  potential  temperature  and  water  vapor 
balances  to  yield 

ae  _  ae  sTo1  i  e  dRk  1C, 

dt  J  dxj  dxj  p  cp  T  dxk 

and 


dq  _  dq  du 'q' 

—  +  U .  — 2-  = - - — 

dt  1  dXj  dxj 


(3.16) 


3.3.  Equations  for  Turbulent  Kinetic  Energy 

An  expression  for  the  time  rate  of  change  of  the  Reynolds  stress  tensor  may  be 
formed  by  subtracting  (3.12)  from  (3.11)  and  applying  the  hydrostatic  and  Boussinesq 
approximations: 


—  6u'  .  du :  ,  du' 

+  U  , - +  U  , - !-  +  U  , - - 

1  a c,  J  ac ,  J  ac 


1  dp'  T;  .  d2u\  du^du'j 

—■ —  +  =  (?8n  +  v - —  + - - 

p  ac,  Tv  dxjdxj  dxj 


The  relation  for  u'p  can  be  written 


(3.17) 


—  du'P  ,  duP  ,  fo'p 

'  +  —  +  W0  ~ -  +  W  — 

dx  qdx  dx 


1  dp'  T;  _  d2u'p  du'du' 
P  dx  Tv  dx  dx  dx 


(3.18) 


Note  that  the  Coriolis  term  is  not  included  because  scale  analysis  shows  that  turbulent 
velocity  fluctuations  occur  at  time  scales  too  short  to  be  influenced  by  the  Coriolis  force 
(Shaw,  1990).  Multiplying  (3.17)  by  u'p  and  (3.18)  by  u\ ,  averaging  both  equations  and 
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adding  them  yields  the  equation  for  the  time  rate  of  change  of  the  Reynolds  stress  tensor 


\{u'pdp'  +  u]dx ^ 


dtu'u” =  UpU-  &; ~ 


dx,  dx 


j  ) 


+  JL(uX83i+u'pt;8:„p)+vi 


fd2u'u'p  2 du ;  du'p  ^ 


^  dxjdxj  dx j  dXj  j 


(3-19) 


Using  the  notation  in  Shaw  (1990),  the  mean  kinetic  energy  of  turbulence  (TKE)  per 
unit  volume  is  defined  as 


1 


E  =  -pu'iu'i 


Alternately,  meteorologists  generally  use  the  kinematic  form  to  define  the  TKE 


-  E  1  — 
e  =  —  =  —  u,ui 
p  2  ' 


and  the  instantaneous  net  contribution  to  the  TKE 


(3.20) 


(3.21) 


e=~u]u’ 


(3.22) 


The  equation  for  rate  of  change  of  the  turbulence  kinetic  energy  follows  directly  from 
(3.19)  by  setting  p=i  and  making  use  of  the  continuity  equation: 


dt  1  dXj  dXj  Tv  dXj  p  dxj  dx ,  dx. 


(3.23) 


The  first  term  of  the  right  hand  side  of  the  equation  represents  the  rate  of  change  of 
TKE  due  to  advection.  The  second  term  represents  the  generation  of  mechanical 
turbulence  through  shear,  while  the  third  term  represents  buoyancy  production.  The 
fourth  and  fifth  terms  are  transport  terms  which  neither  produce  nor  consume  energy,  but 
act  to  redistribute  TKE  within  the  boundary  layer;  the  first  transport  term  represents  the 
divergence  of  turbulent  kinetic  energy  flux,  and  the  second  is  due  to  the  divergence  of 
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pressure-velocity  covariance.  The  last  term  is  always  negative  and  represents  the 
dissipation  a  of  turbulent  kinetic  energy  at  the  smallest  scales  due  to  viscous  forces. 

The  energy  budget  simplifies  greatly  under  horizontally  homogeneous  conditions  for 
which  mean  quantities  depend  only  on  the  vertical  coordinate.  Furthermore,  since  we 
expect  the  smallest  scales  of  motion  to  have  an  isotropic  structure,  the  expression  for 
dissipation  can  also  be  simplified.  Replacing  w,  and  xt  by  (u,v,w)  and  (x,y^), 
respectively,  in  equation  (3.23),  the  TKE  budget  equation  for  stationary,  horizontally 
homogeneous  conditions  with  zero  subduction  (w  =  0)  becomes 


—  uw 


7  du 


dz 

This  is  the  form 


7v^+i77v4^-— 

dz  Tv  p  dz  dz 


of  the  TKE  budget  that  will  be  used 


-e  =  0  (3.24) 

in  the  following  discussions. 
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4.  Methodology 


There  are  four  standard  methods  for  determining  the  momentum  flux  from 
atmospheric  information: 

(i)  direct  covariance  (eddy  correlation); 

(ii)  mean  profiles; 

(iii)  bulk  aerodynamic;  and, 

(iv)  direct  /  inertial  dissipation. 

4. 1 .  Direct  Techniques  to  Determine  the  Momentum  Flux 

4.1.1.  Direct  Covariance  (Eddy  Correlation)  Method 

The  direct  covariance  method  is  the  standard  upon  which  all  other  methods  of 
estimating  the  momentum  flux  are  based.  Whereas  the  other  three  methods  rely  on 
similarity  theory,  the  covariance  method  is  a  direct  computation  of  the  Reynolds  stress 
and  yields  an  unbiased  estimate  of  the  flux  (Fairall  et  al.,  1990).  Revisiting  equation 
(3.14),  the  wind  stress  can  be  determined  directly  by  measuring  the  turbulent  components 
of  the  streamwise  and  vertical  wind  velocity: 

x  =  -p[/  mV  +  ;'vv]  (4. 1 ) 

In  practice,  the  ensemble  average  is  approximated  by  averaging  simultaneous 
measurements  of  the  turbulent  components  of  wind  velocity  and  cross  correlating  them 
over  a  finite  time  interval  (Frederickson  et  al.,  1997).  The  assumption  that  allows  us  to 
relate  time  averaged  fluxes  to  ensemble  averaged  fluxes  is  referred  to  as  ergodicity 
(Panofsky  and  Dutton,  1984).  For  stationary  conditions,  as  the  length  of  the  averaging 
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interval  increases,  the  accuracy  of  this  approximation  generally  increases,  although  there 
is  a  tradeoff  between  the  choice  of  averaging  period  and  nonstationarity  of  the  data 
(Mahrteta/.,  1996). 

The  primary  disadvantage  of  the  direct  covariance  method  is  that  comparative  studies 
(e.g.,  Fairall  et  al.,  1990;  Edson  et  al.,  1991;  Edson  and  Fairall,  1998)  indicate  that  the 
direct  covariance  method  is  highly  susceptible  to  flow  distortion.  Additionally,  unless  the 
measurements  are  made  from  a  stable  platform,  such  as  a  fixed  mast  or  tower,  the  data 
must  be  corrected  for  platform  motion  contamination.  Clearly,  this  limits  the 
applicability  of  this  method  from  ships  and  buoys. 

4.1.2.  Motion-Corrected  Covariance  Method 

However,  a  great  deal  of  effort  has  been  made  in  recent  years  to  compute  direct 
covariance  fluxes  from  moving  platforms  to  improve  flux  parameterizations  over  the 
open  ocean  (e.g.,  Fujitani,  1981,  1985;  Tsukamoto,  1990;  Dugan  et  al.,  1991;  Anctil 
et  al,  1994;  Fairall  et  al.,  1997;  Edson  et  al.,  1998).  If  one  wishes  to  apply  the  direct 
covariance  method  on  a  moving  platform,  then  the  orientation  and  motion  of  the  platform 
must  be  accounted  for  before  computing  the  correlations  (Fairall  et  al,  1990).  Surface¬ 
following  buoys  experience  a  variety  of  motions  including  instantaneous  tilt  of  the 
sensors  due  to  platform  movement,  angular  velocities  about  the  platform’s  local 
coordinate  system,  and  translational  velocities  with  respect  to  a  fixed  frame  of  reference 
(Fujitani,  1981;  Hare  et  al.,  1992).  In  vector  notation,  the  motion  corrected  wind  speed 
Vlrue  can  be  written  as 

^=T(ni,+£2Jx«)+f„  (4.2) 
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where  T  is  the  rotation  matrix  which  transforms  measurements  taken  in  the  buoy 
reference  frame  to  earth  coordinates;  Vobs  and  fiA  are  the  measured  wind  speeds  and 

angular  rates  in  the  buoy  frame  of  reference,  respectively;  R  is  the  position  vector  of  the 
wind  sensor  with  respect  to  the  motion  package;  and,  Vmot  represents  the  translational 
velocities  in  the  fixed  frame  (Anctil  et  al.,  1994;  Edson  et  al.,  1998).  The  product 

x  TR  represents  the  angular  velocity  in  the  fixed  frame. 

Removal  of  this  motion  from  our  wind  sensors  requires  us  to  measure  the  three 
dimensional  translational  motion  of  a  reference  point  on  the  buoy,  the  three  angular 
rotations  about  that  point,  and  the  Euler  angles  required  for  the  transformation  matrix. 
During  the  CMO  experiment,  three  Systron  Donner  quartz  sensors  recorded  the  pitch, 
roll,  and  yaw  rates  (0,0,  and  'P ,  respectively),  and  a  three-axis  Columbia  linear  servo 
accelerometer  provided  the  linear  accelerations  of  the  buoy  ( x  ,y ,  and  z ,  respectively)  . 
All  of  the  sensors  were  strapped  down  in  the  hull  of  the  discus  buoy.  One  of  the  main 
advantages  of  a  strap-down  system  is  that  it  readily  computes  the  angular  rates  in  the 
buoy  frame  of  reference: 

Qb  =  [®i  +  0j  +  'f'k]  (4.3) 

where  i ,  ] ,  and  k  represent  the  unit  vectors  in  a  right-hand  coordinate  system  pointing 
from  the  buoy  fin  to  front,  starboard  to  port,  and  vertically  up,  respectively. 

The  three  Euler  angles  which  specify  the  rotation  of  the  buoy  are  the  pitch  O,  roll  0, 
and  yaw  *F.  Using  the  approach  described  by  Edson  et  al.  (1998),  the  Euler  angles  are 
approximated  using  complimentary  filtering,  which  combines  the  high-passed  filtered 
integrated  angular  rates  with  the  low-passed  linear  accelerations  as 
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(4.4) 


® 


1  y  is  _ 

■  + - ® 


x  + 1  g  x  s  + 1 


..  lx  TS  _ 

0  « - 1 - © 

X  + 1  g  X 5  +  1 


(4.5) 


and 


1 

x  +  1 


s/ow 


+-^-vp 

xs  +  l 


(4.6) 


where  s  represents  the  differentiation  operator  such  that  ®  =  5® ,  x  is  a  time  constant, 
^/ow  represents  the  compass  output,  and  a  small  angle  approximation  is  made  in  (4.4) 
and  (4.5). 

Given  the  Eulerian  angles,  the  rotation  matrix  T  which  transforms  measurements 
taken  in  the  buoy  reference  frame  to  earth  coordinates  can  be  expressed: 


T  = 


cos  *F  cos© 
sin  4*  cos  0 


-  sin© 


-  sin  *F  cos®  +  cos  *F  sin  ©sin®  sin 'F  sin®  +  cos 'F  sin  ©cos® 
cos  'F  cos  ®  +  sin  'F  sin  0  sin  ®  sin  0  cos  ®  sin  'F  -  sin  ®  cos  *F 

cos  0  sin  ®  cos  0  cos  ® 


(4.7) 


Again,  the  sign  convention  corresponds  to  a  right-handed  coordinate  system  with1? 
positive  for  the  ship’s  bow  yawed  counter-clockwise  from  north,  ®  positive  for  the  port 
side  rolled  up,  and  ©  positive  for  the  bow  pitched  down  (Edson  et  al.,  1998). 

The  accelerations  are  then  rotated  into  the  vertical  using  the  rotation  matrix  T  before 
adding  the  gravity  vector,  g  =  (0,0, -g) ,  to  the  rotated  accelerations.  The  resulting 
accelerations  are  then  high-pass  filtered  to  remove  any  drift  before  integrating.  This 
process  can  be  expressed  as: 

W.«  =  §Hp(n(0  +  g)}*  (4.8) 
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where  Hp  represents  a  high-pass  filter  operator.  The  strapped-down  system  thereby 
provides  all  the  terms  required  by  (4.2)  to  compute  the  true  wind  velocity.  These 
corrected  velocities  are  then  used  to  compute  the  wind  stress. 

4.2.  Indirect  Techniques  to  Determine  the  Momentum  Flux 

Monin  and  Obukhov  (1954)  postulated  that  the  turbulent  fluxes  of  atmospheric 
momentum  and  heat  energy  in  a  stationary,  horizontally  homogeneous  surface  layer  is 

determined  by  the  height  above  the  surface  z,  the  buoyancy  parameter  g/Qv  ,  the  friction 
velocity  u. ,  and  the  surface  buoyancy  flux  w'Q'V0  where  %v  is  the  virtual  potential 

temperature  and  accounts  for  density  changes  due  to  humidity  fluctuations  and  g  is  the 
acceleration  of  gravity  (Lumley  and  Panofsky,  1963;  Wyngaard,  1973).  The  friction 
velocity  u,  is  defined  in  the  manner  of  equation  (4.1)  from  the  surface  stress: 

=  [(-iiV)2+(-vV)2f  (4.9) 

Largely  on  the  basis  of  dimensional  arguments,  Monin-Obukhov  (MO)  similarity 
theory  predicts  that  the  various  terms  in  the  TKE  budget  [equation  (3.24)]  are  universal 
functions  of  the  MO  stability  parameter  C,  after  normalization  by  the  appropriate  scaling 

parameter,  — j : 

tit 

u, 

where  k  is  the  dimensionless  von  Karman  constant,  generally  assumed  to  be  0.4 
(Hdgstrom,  1987;  Edson  et  al,  1991;  Fairall  et  al.,  1996a);  <j»E  is  the  dimensionless 
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dissipation  function;  §tp  and  <J>/e  are  the  dimensionless  transport  terms;  and,  <|>M  is  the 
dimensionless  shear  defined  as 


kz  (du  V  ( dv  V 

<MC)  =  —  -r-  +  t- 

ut  \oz)  \oz  J 


(4.11) 


The  stability  parameter  ^  represents  the  ratio  of  the  height  above  the  surface  z  to  the 

q  u  3  2 

Monin-Obukhov  length  L  defined  as  L  =  — — —  *  ,  such  that  C  =  — .  In  unstable 

gKW0,o  s  L 

conditions  the  Monin-Obukhov  length  represents  the  height  where  mechanical  and 
buoyant  production  of  turbulence  are  approximately  equal,  i.e. ,  at  this  height,  C,  =  -1 .  At 
heights  below  this  level,  shear  production  dominates;  while  above  this  level,  buoyancy 
dominates.  In  stable  conditions  the  stratification  consumes  some  of  the  kinetic  energy 
such  that  ^  represents  the  ratio  of  buoyant  consumption  to  shear  production. 

4.2. 1 .  Mean  Profile  Method 

Generally,  the  wind  direction  is  assumed  to  be  constant  in  the  surface  layer  such  that 
the  coordinate  system  is  rotated  so  that  v  =  0 .  In  this  coordinate  system  the  wind  shear 
becomes 


OZ  KZ 


(4.12) 


Upon  integrating  the  above  equation,  one  obtains  the  logarithmic  wind  profile  for  neutral 
conditions  upon  which  the  profile  method  is  based: 


w.  .  z  z 

m(z)  =  m(z0)  +  —  In - V|/M  - 


(4.13) 
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where  z0  is  the  aerodynamic  roughness  length  and  reflects  the  roughness  of  the  surface. 


and  v|/M  is  the  integrated  form  of  <|>M  .  Using  equation  (4.13),  it  is  relatively  easy  to 
obtain  w,  and  z0  from  a  given  wind  profile:  plotting  u  as  a  function  of  ln(z)  -y|/m(Q 

results  in  a  straight  line  with  —  as  the  slope  and  —  ln(z0)  as  the  intercept. 

K  K 

Therefore,  the  implementation  of  the  profile  method  requires  the  deployment  of 
multiple  sensors  at  varying  heights.  This  type  of  deployment  in  the  marine  surface  layer 
creates  a  number  of  problems.  In  addition  to  the  cost  associated  with  employing  an  array 
of  sensors,  excellent  relative  calibration  among  sensors  is  required  due  to  the  small  shear 
found  over  the  ocean.  Additionally,  calibration  drift  between  the  individual  sensors  due 
to  mechanical  degradation  may  result  in  large  uncertainties  among  measurements. 
Another  physical  limitation  of  the  profile  method  is  the  requirement  for  the  lowest  sensor 
to  be  at  least  two  or  three  wave  heights  above  the  interface  since  the  presence  of 
undulating  long  waves  causes  deviations  to  the  logarithmic  profile  (Geemaert,  1990). 
This  requirement  further  limits  the  placement  of  the  sensors  since  they  upper  height  is 
already  bounded  by  the  height  of  the  surface  layer,  the  region  within  which  the  fluxes 
vary  by  no  more  than  10%  of  their  surface  values,  a  requirement  of  MO  similarity  theory. 
Finally,  and  probably  most  importantly,  the  physical  set-up  of  the  array  renders  the 
sensors  extremely  sensitive  to  flow  distortion  and  motion  contamination.  Therefore,  this 
method  has  rarely  been  implemented  over  the  ocean. 
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4.2.2.  Bulk  Aerodynamic  Method 

The  bulk  method  is  the  simplest  and  most  widely  used  of  the  four  techniques.  The 
bulk  method  is  an  application  of  the  profile  method  for  which  the  structure  of  the 
atmospheric  boundary  layer  is  estimated  by  using  measurements  at  one  level  and 
assuming  a  logarithmic  wind  profile.  The  standard  bulk  expression  for  the  stress 
components  is 

T0;  =p«.2  =paCD(M/(2)-M/(z0)]M,(2)-M;(z0)|  (4.14) 

where  CD  is  the  bulk  transfer  coefficient  for  stress;  u,  ( z )  is  the  average  of  one  of  the 
horizontal  wind  components  relative  to  the  fixed  Earth  measured  at  some  reference  height 
z\  and  m,  (z0)  is  the  surface  current  (Fairall  et  al.,  1996a).  The  velocity  difference  is, 


therefore,  the  velocity  relative  to  the  ocean  surface. 

Equation  (4.14)  can  then  be  used  to  define  the  drag  coefficient  as 


Q(z)  = 


n2 


K 


[log(z/zo)-tMOj 


(4.15) 


The  neutral  drag  coefficient  is  related  to  the  roughness  length  by 


CdN  (z)  — 


log  (z/z0). 


(4.16) 


such  that 


X 


1  + 


Vm(Q 


(C„f 


(4.17) 


K 

where  the  subscript  N  denotes  the  value  in  neutral  conditions;  and,  \\iM  (^)  is  the  MO 
similarity  profile  function.  The  following  expressions  of  Businger  et  al.  (1971)  and  Dyer 
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(1974)  are  the  most  commonly  used: 


<Kj<(?)=(1-1«S)*  ?<»  (4.18a) 

VvfeMi+s?)  c>0  (4.18b) 

Using  the  approach  described  by  Smith  (1988),  the  velocity  aerodynamic  roughness 
length  z0  can  be  expressed  as 


z,=zs+zr  (4.19) 

where  zs  and  zr  represent  the  aerodynamic  roughness  of  smooth  and  rough  flows, 

respectively.  The  roughness  length  for  a  smooth  surface  depends  on  the  viscosity  of  air 
v  and  the  friction  velocity  (Businger,  1973): 


(4.20) 


The  determination  of  zr  is  based  on  the  relationship  proposed  by  Chamock  (1955). 

Using  dimensional  analysis,  Chamock  assumed  that  the  aerodynamic  roughness  length  of 
the  surface  for  well-developed  seas  is  proportional  to  the  wind  stress  which  generates  and 
supports  the  surface  waves  and  gravity  which  acts  as  the  restoring  force: 


8 


where  ac  is  the  Chamock  “constant”  for  which  values  between  0.010  and  0.035  can  be 

found  in  the  literature  (Garratt,  1992).  Representing  the  ratio  of  inertial  force  to 
gravitational  force,  the  Chamock  constant  is  analogous  to  a  Froude  number. 

The  method  and  algorithm  for  computing  the  bulk  estimates  are  described  in  detail  in 
Fairall  et  al.  (1996a).  In  summary,  the  bulk  estimates  are  calculated  by  using  the 
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roughness  length  described  above  with  a  Chamock  constant  of  0.01 1  and  with  the 
following  corrections: 

(i)  In  highly  unstable  conditions,  £  =  f  «  -1 ,  the  profiles  of  temperature, 
moisture, 

and  momentum  are  altered  so  that  they  obey  the  C,  ^  asymptotic  convective  limit 
dependence  on  stability  (Wyngaard,  1973). 

(ii)  A  gustiness  velocity  wg  is  used  in  the  computation  of  the  mean  wind 

speed,  p  is  an  empirical  constant  of  the  order  of  1 .0,  and  the  convective  velocity  W.  is 
computed  as  follows: 


W ?  =- 


Qh 

pac 


0.61  T- 


pa 


P  A 


(4.22) 


where  QH  and  QE  are  the  sensible  and  latent  heat  fluxes,  respectively;  pa  and  are 


the  density  and  isobaric  specific  heat  of  air;  Le  is  the  latent  heat  of  vaporization;  T  is  the 
potential  temperature  at  the  surface;  and  zt  is  the  depth  of  the  convective  boundary  layer 

(approximately  600  m).  This  factor  accounts  for  the  appreciable  momentum  flux 
observed  in  highly  convective  conditions  where  the  mean  wind  shear  goes  to  zero  but  its 
variance  (i.e.,  gustiness)  does  not. 

(iii)  The  “Webb”  correction  (Webb  et  al.,  1980),  which  accounts  for  variations  in  air 
density  due  to  the  presence  of  a  heat  flux,  is  applied. 

(iv)  Cool-skin  correction  as  described  by  Fairall  et  al.  (1996b)  for  nonradiative  sea 
surface  temperature  measurements  is  incorporated.  Note  that  since  the  temperature 
sensor  was  at  1 .5  m,  no  warm  layer  correction  was  applied  due  to  the  shallow  depth. 
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One  potential  source  of  error  in  the  bulk  aerodynamic  method  is  in  the 
characterization  of  the  roughness  length  by  equations  (4.19)  through  (4.21).  Chamock’s 
expression  for  the  roughness  of  the  sea  surface  has  been  shown  to  be  generally  valid  for 
moderate  to  high  wind  speeds  at  unlimited  fetch  and  duration  (Geemaert  et  al.,  1986). 
However,  upon  examining  equation  (4.21),  one  sees  that  no  characteristics  of  the  wave 
field  are  truly  present.  This  implicitly  assumes  that  the  wind  and  wave  fields  are  in  local 
equilibrium;  consequently,  Chamock’s  relation  will  only  give  consistent  results  for  a 
fully-developed  sea  (Donelan,  1990;  Maat  et  al.,  1991;  Dobson  et  al.,  1994).  Sudden 
changes  in  the  wind,  such  as  that  due  to  frontal  passages,  result  in  the  wind  field  being 
temporarily  out  of  equilibrium  with  the  existing  wave  field  (Smith,  1988).  The 
momentum  transfer  between  the  wind  and  wave  fields  as  the  waves  grow  or  decay  in 
response  to  the  forcing  is  not  considered  in  Chamock’s  expression.  Additionally,  the 
influence  of  varying  fetch  and  possible  interactions  with  the  bottom  in  shallow  water  are 
not  included.  The  neglect  of  these  processes  probably  explains  the  wide  range  of  values 
for  the  Chamock’s  “constant”  found  from  different  experiments.  Finally,  implementation 
of  the  bulk  aerodynamic  method  using  buoy  data  often  places  the  measurements  where 
the  waves  cause  changes  to  the  log  profile.  All  of  these  effects  result  in  changes  to  the 
drag  coefficient  which  are  sea-state  dependent. 

4.2.3.  Direct  and  Inertial  Dissipation  Methods 

The  direct  dissipation  and  the  inertial  dissipation  methods  differ  only  in  the  manner 
in  which  the  dissipation  rate  for  turbulent  kinetic  energy  is  computed.  Direct 
measurements  of  dissipation  rates  are  determined  by  integrating  over  the  spectrum  of  the 
velocity  derivatives  (Fairall  and  Larsen,  1986): 
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Since  the  fluctuations  have  a  spatial  resolution  that  approaches  the  Kolmogorov 
microscale  (typically  1  mm  in  atmospheric  flows),  hot-wire  anemometers  must  be  used  to 
obtain  direct  measurements  (Edson  and  Fairall,  1998).  While  this  method  has  been  used 
overland,  its  applicability  over  the  ocean  is  limited  due  to  contamination  and  even 
destruction  of  the  wires  by  sea  spray  (Edson  and  Fairall,  1998).  Therefore,  the  direct 
dissipation  method  is  rarely  employed  over  the  ocean. 

Instead,  the  inertial  subrange  method,  which  requires  frequency  responses  only  up  to 
the  order  of  10-100  Hz,  is  more  commonly  employed  (Smith  et  al.,  1996;  Fairall  et  al., 
1997).  In  the  inertial  subrange  of  local  isotropic  turbulence,  the  one-dimensional  velocity 
variance  spectrum  Ou  can  be  expressed  as  a  function  of  wavenumber  magnitude  as 

<bu(k)  =  auz%k~SA  (4.24) 

where  au  is  the  one-dimensional  Kolmogorov  constant,  taken  to  be  0.52  (Hogstrom, 


1996),  and  k  is  the  wavenumber. 

Invoking  Taylor’s  frozen  turbulence  hypothesis,  k  =  ^==-  s  = ,  the  velocity 


variance  spectra  can  be  related  to  the  frequency  spectra  Su  as 


M>„(*)  = 


A(/) 


+  att2 ) 


=  aus 
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'2g_V/ 

kU.j 


(4.25) 


where /is  frequency,  Su  is  the  spectrum  of  the  streamwise  velocity  Us  as  a  function  of 
the  angular  frequency  to ,  and  Tuu  and  (att]  +  att2 )  are  correction  factors: 


(i)  Tuu  is  a  correction  for  inaccuracies  in  using  Taylor’s  hypothesis  to  estimate  the 


magnitude  of  wavenumber  spectra  in  the  inertial  subrange  and  has  the  form  given  by 
Wyngaard  and  Clifford  (1977): 


=  1  ik)2  ,2(gv+^) 

9  U }  3  U } 


(4.26) 


2  2  2 

where  au  ,  av  ,  aw  are  variances  (uu,  w  and  ww)  computed  from  the  same  time  series 
that  produced  the  spectral  estimate  (Hill,  1996). 

(ii)  The  terms  attx  and  att2  correct  for  excess  damping  of  higher  frequencies  caused 
by  the  eight-fold  pulse-averaging  of  the  sonic  anemometer  (Henjes  et  al.,  1998)  and  has 
the  forms 
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(4.27b) 


where  /0  is  the  Nyquist  frequency,  20.833  Hz. 

Note  that  spatial  averaging  also  reduces  the  spectral  estimates  at  high  frequencies,  in 
the  manner  described  by  Kaimal  et  al.  (1968): 
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(4.28) 


where  /  is  the  separation  distance  between  the  sonic  heads.  However,  the  effect  of  line¬ 
averaging  is  negligible  for  kl  <  1 .  Since  we  have  limited  the  frequency  range  such  that 

/  <  — ,  this  correction  is  not  needed. 

J  2nl 

When  using  the  inertial  dissipation  method,  it  is  generally  assumed  that  the  transport 
terms  cancel  (Hicks  and  Dyer,  1972;  Large  and  Pond,  1981;  Dyer  and  Hicks,  1982;  Large 
and  Businger,  1988)  so  that  a  balance  exists  between  the  viscous  dissipation  and 
production  (buoyant  and  shear)  of  TKE.  Under  this  assumption,  equation  (4.10)  reduces 
to 

“ (4.29) 

w. 

A  potential  source  of  error  in  using  the  inertial  dissipation  technique  arises  from  the 
assumption  that  there  is  a  balance  between  TKE  production  and  dissipation.  Recent 
investigations  (Vogel  and  Frenzen,  1992;  Frenzen  and  Vogel,  1992;  Thiermann  and 
Grassl,  1992;  Yelland  and  Taylor,  1996;  Edson  et  al.,  1997;  Edson  and  Fairall,  1998) 
suggest  that  there  is  a  local  imbalance  between  the  pressure  and  energy  transport  terms 
such  that  the  magnitude  of  the  energy  transport  term  <|>te  (C,)  is  greater  than  the  magnitude 

of  the  pressure  transport  term  <|>,p(c),  and  TKE  production  exceeds  dissipation  by  15- 

20%.  Over  the  ocean,  it  has  been  hypothesized  by  Edson  et  al.  (1997)  that  this 
dissipation  deficit  is  due  to  the  flux  of  energy  going  into  ocean  waves. 
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The  flux  of  kinetic  energy,  F(z),  into  a  layer  of  air  over  a  horizontally  homogeneous 
surface  is  given  by 

F(z)  =  u'w'u(z)  +  v'w'v(z)  +  w'e'  +  —  xv'p'  (4.30) 

P 

where  the  first  two  terms  on  the  right-hand  side  of  the  equation  represent  the  flux  of  mean 
flow  kinetic  energy  and  the  last  two  terms  represent  the  rate  of  diffusion  of  kinetic  energy 
(Edson  et  al.,  1997).  Over  land,  it  is  generally  assumed  that  there  is  no  energy  flux 
through  the  ground  so  that  the  flux  entering  the  layer  at  height  h  can  be  related  to  the  total 
rate  of  dissipation  within  the  layer  by 

h  h _ 

Je  dz  =  -F(h)  +  —  Jw'  T'dz  (4.31) 

0  0 

where  the  last  term  represents  the  generation  of  kinetic  energy  due  to  buoyancy  flux 
(Edson  et  al.,  1997).  For  neutral  conditions,  equation  (4.3 1)  states  that  the  flux  of  kinetic 
energy  into  the  layer  equals  the  dissipation  within  that  layer;  i.  e. ,  production  balances 
dissipation.  However,  over  the  ocean,  where  the  exchange  of  momentum  from  the 
atmosphere  to  the  ocean  results  in  wind-driven  waves  and  currents,  the  surface  energy 
flux  F{ 0)  is  not  negligible  and  the  additional  forcing  mechanism  of  the  surface  waves 
must  be  included. 

The  effects  of  surface  waves  are  included  in  the  KE  energy  budget  equation  (3.24) 
and  energy  flux  equation  (4.30)  by  decomposing  the  velocity  into  mean,  turbulent,  and 
wave-induced  components 

u(t)  =  u  +u'(t)  +  u(t)  (4.32) 
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where  the  tilde  represents  the  wave-induced  component.  The  wave-induced  component  is 
defined  with  the  use  of  the  phase-averaging  operator  as  described  by  Finnegan  et  al. 
(1984): 


(w(O)p  =w+M(0  =  lim 

Vi->  0 


Nl 

£«( t  +  mo ) 

n= 1 


«(0  =  («(0)F  -« 


(4.33) 

(4.34) 


where  xa  is  the  period  of  the  organized  oscillation  and  NL  represents  the  number  of 

consecutive  wavelengths  averaged  together  to  get  one  phase-averaged  wave.  This  type  of 
averaging  can  be  used  to  rewrite  equation  (3.24)  to  include  the  wave-induced 
contributions.  An  approximation  of  the  steady  state  solution  of  the  TKE  energy  budget  is 
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=  0 


dz  dz 

where  rn  and  rtj  can  be  viewed  as  the  wavelike  fluctuations  in  the  background  Reynolds 


stress  due  to  the  presence  of  waves 

={u\u'J)F-uyj  (4.36) 

such  that  the  sum  of  the  last  two  terms  in  brackets  represents  the  interaction  between  the 
wave  field  and  the  wavelike  part  of  the  turbulent  stresses  and  constitutes  a  loss  of  kinetic 
energy  (Finnigan  et  al.,  1984).  In  this  form,  the  shear  production  P  now  includes  terms 
representing  the  wave-induced  momentum  flux, 

p  =  +  +  (4.37) 

dz  dz 
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The  production  of  kinetic  energy  is  therefore  the  result  of  a  forcing  mechanism  that  is  a 
function  of  sea  state. 

Similarly,  the  transport  of  energy  from  wind  to  waves  at  the  ocean  surface  is  actually 
the  lower  boundary  condition  of  the  pressure  transport  term,  i.  e. , 

f(0)  =  “H7>L*“A,  (4.38) 

p  p  ot 

where  %  is  the  surface  elevation  and  the  time  rate  of  change  is  known  as  the  piston 


velocity.  This  term  must  be  included  in  equation  (4.31)  in  over  ocean  studies  such  that 


(4.39) 


where  the  fluxes  now  represents  both  turbulent  and  wave-induced  fluctuations.  Unlike 
the  solid  boundary  for  which  the  flux  through  its  surface  is  negligible,  there  is  a  transfer 
of  energy  from  the  atmosphere  to  the  ocean,  and  the  total  amount  of  dissipation  within  the 
wave  boundary  layer  is  expected  to  be  less  than  required  to  balance  the  energy  flux  into 
the  layer;  a  phenomenon  supported  by  recent  investigations  (Yelland  and  Taylor,  1996; 
Edson  et  al. ,  1997;  Edson  and  Fairall,  1998).  Thus,  in  the  wave  boundary  layer,  the 
transport  of  energy  and  momentum  to  the  ocean  from  the  atmosphere  is  not  adequately 
described  using  Monin-Obukhov  similarity  theory,  and  additional  scaling  parameters  are 
required  for  similarity. 
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5 .  Analysis  and  Results 


The  objective  of  this  analysis  is  to  improve  estimates  of  the  momentum  flux  using 
modified  forms  of  the  bulk  aerodynamic  and  inertial  dissipation  methods.  The 
modifications  to  the  bulk  aerodynamic  method  involve  the  derivation  of  a  roughness 
length  parameterization  which  is  sensitive  to  sea  state.  These  modifications  are  drawn 
from  comparisons  with  direct  covariance  and  inertial  dissipation  flux  estimates.  In  the 
case  of  the  inertial  dissipation  momentum  fluxes,  additional  parameterizations  are  applied 
to  account  for  the  effects  of  surface  waves  and  the  non-zero  energy  flux  from  the 
atmosphere  to  the  ocean. 

5.1.  Data  Selection 

The  data  sets  obtained  from  both  the  VAWR  sensor  and  sonic  anemometer  span  the 
entire  eleven  months  from  August  1996  to  June  1997  and  represent  a  total  of  15,262  data 
points  (2  measurements  per  hour  for  approximately  3 1 8  days).  Of  these  measurements  a 
total  of  2,303  points,  or  15%,  were  removed  prior  to  processing  for  one  or  more  of  the 
following  reasons: 

(i)  Instrument  failure  (VAWR  and/or  sonic  anemometer). 

(ii)  The  velocity  spectrum  for  inertial  dissipation  calculations  did  not  meet  the  -  f 
slope  criteria  [see  equation  (4.24)]  and  did  not  allow  the  calculation  of  dissipation 
estimates.  Possible  reasons  include  rain  and  spray  contamination  and  obstructions  to  the 
path  (e.g.,  due  to  a  seagull  sitting  on  the  sensor). 

(iii)  Flow  distortion,  which  arises  when  the  relative  wind  direction  is  greater  than 
±  60°  from  head  on. 


39 


The  wave  height  data  set  represents  241  days  of  the  eleven-month  deployment,  but 
only  1  measurement  per  hour  was  taken  for  a  total  of  5,778  points.  In  order  to  insure  that 
oceanographic  conditions  resulted  from  locally-generated  seas  and  not  from  swell 
propagating  into  the  test  site,  8%  of  the  wave  data  were  discarded  prior  to  processing 
because  two  distinct  peaks  were  present. 

For  motion-corrected  covariance  calculations,  the  data  set  was  further  limited  by  the 
amount  of  data  collected  by  the  motion-correction  package.  This  package  collected  30 
minutes  of  1 0  Hz  data  at  noon  and  midnight  for  3 1 7  days  which  produced  a  data  set 
totaling  634  direct  covariance  flux  estimates. 

5.2.  Traditional  Approach 

Using  the  procedures  outlined  in  Chapter  4,  the  friction  velocity  u. ,  as  defined  by 
equation  (4.9),  was  calculated  using  the  bulk  aerodynamic,  inertial  dissipation,  and  direct 
covariance  methods.  Figures  5.1  through  5.3  are  comparisons  of  the  u,  calculations.  The 
inertial  dissipation  and  bulk  aerodynamic  methods  provide  remarkably  similar  values  for 
u,  <  0.6  m/s  (Figure  5.1).  At  higher  wind  speeds,  the  inertial  dissipation  method  predicts 
progressively  higher  values  of  u,  than  the  bulk  aerodynamic  method,  although  the  values 
generally  exhibit  less  than  a  12%  difference.  Compared  to  friction  velocity  values 
computed  by  the  motion-corrected  covariance  method  (Figures  5.2  and  5.3),  both  the  bulk 
aerodynamic  and  inertial  dissipation  methods  provide  u.  estimates  which  are  consistently 
7-12%  lower  than  expected  for  w,  >  0.2  m/s.  This  would  suggest  that  although  the  two 
methods  provide  quite  similar  results,  these  values  are  underestimating  the  actual 
momentum  flux  over  the  ocean. 
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Bulk  Aerodynamic  u.  (m/s) 


Friction  Velocity  u*  Comparison 


Figure  5.1  A  comparison  of  the  friction  velocities,  u, ,  estimated  using  the  bulk 

aerodynamic  versus  inertial  dissipation  methods  and  bin-averaged  in  increments  of 
0.02  m/s.  The  error  bars  denote  the  standard  error  (standard  deviation  divided  by  the 
number  of  points). 


Original  Bulk  Aerodynamic  Friction  Velocity  u*  Comparison 


Figure  5.2  Friction  velocity,  ut ,  comparisons  of  the  bulk  aerodynamic  versus  motion- 
corrected  covariance  method.  Values  of  u,  are  bin-averaged  in  increments  of 
0.02  m/s. 
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Original  Inertial  Dissipation  Friction  Velocity  u*  Comparison 


Figure  5.3  Friction  velocity,  ut ,  comparisons  of  the  inertial  dissipation  versus  motion- 
corrected  covariance  method.  Values  of  u„  are  bin-averaged  in  increments  of 
0.02  m/s. 
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The  neutral  drag  coefficient  CDN  described  by  equation  (4.17),  was  calculated  by  the 


bulk  aerodynamic,  inertial  dissipation,  and  direct  covariance  methods  for  the  following 
four  wave-age  bins: 


u. 


C„ 

(ii)  24  <  — <  28 ; 

w* 

C 

(iii)  28  <  —  <  32 ;  and, 

«. 

(iv)  32  <  \ 

u. 


where  Cp  is  the  phase  speed  of  the  dominant  seas,  calculated  using  equation  (2.7).  The 


first  wave-age  bin  corresponds  to  young  (developing)  seas,  the  second  and  third  wave-age 
bins  represent  seas  for  which  the  waves  receive  no  additional  energy  input  from  the  wind 
and  are  considered  folly-developed,  and  the  last  wave-age  bin  corresponds  to  waves  that 


are  decaying.  [An  alternate  expression  for  the  wave  age  is 


— — ,  where  young  waves 
U\0N 


correspond  to  0.22  <  — —  <  0.29  and  folly-developed  waves  are  represented  by 

UX0N 

c 

0.5  <  — —  <1.25  (Donelan,  1990).]  Interestingly,  less  than  10%  of  the  data  represent 
Ul0N 

conditions  in  which  the  seas  were  growing.  The  absence  of  young  waves  from  the  data 
set  most  likely  indicates  a  limit  in  the  techniques  used  to  detect  young  waves.  A  review 
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of  the  meteorological  record  reveals  that  the  CMO  test  site  experienced  numerous  frontal 
passages  and  storm  systems  (including  two  tropical  storms)  with  wind  conditions  that 
should  be  associated  with  new  wave  development.  One  limitation  is  the  inability  of  the 
waverider  buoy  to  sense  waves  with  wavelengths  less  than  twice  its  diameter,  about  2  m. 
Thus,  the  wind  stress  can  only  be  related  to  the  shortest- wavelength  waves  (4  m)  sensed 
by  the  buoy.  Another  possibility  is  that  new  seas  are  masked  by  residual  energy  from 
large  swell  peaks  which  exist  near  the  sea-swell  transition  frequency,  but  this  hypothesis 
has  not  been  explored.  Dobson  et  al.  (1994)  reported  the  same  lack  of  young  seas  in  their 
data  set,  and  both  these  results  indicate  the  need  for  improved  methods  to  detect  and 
categorize  seas.  Nevertheless,  while  we  are  unable  to  distinguish  the  youngest  waves, 
categorization  of  the  seas  into  the  four  wave-age  bins  defined  above  was  possible  and 
provided  invaluable  insight  into  the  relationship  between  the  sea  surface  wind  stress  and 
wave  state. 

The  bulk  aerodynamic  neutral  drag  coefficients,  CDN ,  are  approximately  5-15% 

lower  than  the  values  measured  using  the  direct  covariance  method  for  young  seas  (Figure 
5.4).  For  decaying  seas,  however,  the  bulk  estimates  tend  to  over  predict  the  drag 
coefficient  by  approximately  the  same  error.  As  Figure  5.5  depicts,  the  required 
Chamock’s  “constant”  to  provide  the  best  agreement  for  the  young  seas  is  ac  =  0.0150 ;  a 

value  of  ac  =0.0126  gives  the  best  agreement  for  fully-developed  seas;  while  a  value  of 
ac  =  0.0059  provides  the  best  agreement  for  decaying  seas.  The  Chamock’s  constant 
which  gives  the  best  agreement  when  all  wave-ages  are  considered  is  ac  =0.0116,  which 
is  in  good  agreement  with  the  value  obtained  by  Fairall  et  al.  (1996a)  and  used  in  the 
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Drag  Coefficient  vs  Wind  Speed  Comparison 


C  /u,<=24 
P 


28<Cp/u,<=32 


U10N  if1  ^ 


24<C  /u,<=28 
P 


x  x  Inertial  Dissipation 
+  +  Bulk  Aerodynamic 
o  o  Direct  Covariance 


Figure  5.4  A  comparison  of  the  neutral  drag  coefficient  CDN  estimated  using  the  bulk 
aerodynamic,  inertial  dissipation,  and  motion-corrected  covariance  methods 
(averaged  in  increments  of  UX0N  =  lm/s)  for  four  different  wave-age  ~  bins.  The 
error  bars  for  the  covariance  measurements  denote  the  standard  error. 
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xIO3  C^xlO3 


Bulk  Aerodynamic  &  Direct  Covariance  Drag  Coefficients 


Cp/u.<=24  24<Cp/u.<32 


+  +  Bulk  Aerodynamic 

o  o  Direct  Covariance 

Figure  5.5  A  comparison  of  the  neutral  drag  coefficient,  CDN ,  estimated  using  the  bulk 
aerodynamic  method  and  compared  to  measurements  obtained  by  the  direct 

covariance  method  for  four  different  wave-age  bins.  For  each  wave-age,  the  solid 

2 

oc  u 

curves  represent  the  best  fit  using  Chamock’s  relation,  zr  =  c  *  ,  where  ac 

g 

represents  the  obtained  Chamock’s  “constant”  for  each  wave-age  category.  The  data 
have  been  bin  averaged  in  increments  of  UX0N  =  lm/s,  and  the  error  bars  for  the 
covariance  measurements  denote  the  standard  error. 
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Tropical  Ocean-Global  Atmosphere  Coupled-Ocean  Atmosphere  Response  Experiment 
(TOGA-COARE).  The  results  obtained  here  imply  that  the  drag  coefficient  is  inversely 
correlated  with  wave-age;  i.e.,  the  drag  is  higher  for  young  (developing)  seas  and  lower 
for  fully-developed  seas.  This  result  has  been  theorized  by  Janssen  (1989)  and  supported 
by  experiments  in  the  field  ( e.g .,  Smith,  1980;  Large  and  Pond,  1981;  Geemaert  et  al., 
1986, 1987,  1988;  Davidson  et  al.,  1991;  Maat  et  al.,  1991;  Smith  et  al.,  1992;  Donelan 
et  al.,  1993).  Furthermore,  this  result  supports  the  hypothesis  that  a  roughness  length 
computation  (and,  hence,  drag  coefficient)  which  includes  characteristics  of  the  wave 
field  would  provide  improved  bulk  aerodynamic  momentum  flux  estimates. 

Like  the  bulk  aerodynamic  values,  the  neutral  drag  coefficients  estimated  by  the 
inertial  dissipation  method  do  not  exhibit  a  sea  state  dependence.  For  any  given  wind 
speed,  CDN  is  approximately  the  same  value  regardless  of  the  wave-age.  The  inertial 

dissipation  estimates,  however,  display  a  different  trend  than  the  bulk  values  when 
compared  to  the  direct  covariance  measurements:  the  inertial  dissipation  values 
underestimate  the  actual  drag  for  the  youngest  seas,  but  become  progressively  more 
accurate  as  the  seas  become  fully-developed.  This  implies  that  the  inertial  dissipation 
method  experiences  the  most  error  when  estimating  the  momentum  flux  over  the 
youngest  seas,  but  that  its  performance  improves  with  wave-age.  As  with  the  bulk 
aerodynamic  method,  these  results  suggest  that  incorporating  sea-state  dependence  into 
the  parameterizations  used  in  the  inertial  dissipation  model  would  improve  momentum 
flux  estimates. 
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5.3.  Proposed  Inertial  Dissipation  Method 


The  hypothesis  that  the  dissipation  profiles  are  influenced  by  wave-induced  flow  near 
the  ocean  surface  implies  that  MO-similarity  must  break  down  in  the  wave  boundary 
layer.  This  hypothesis  can  be  investigated  by  comparing  the  vertical  profiles  of  the 
dissipation  rate  with  their  MO-similarity  prediction 

3 


KZ 


(5.1) 


where  C,  =  zjL  and  we  use  the  form  of  the  dimensionless  dissipation  function  given  by 


Edson  and  Fairall  (1998) 

*.(0  =  1 f^-5  ?<»  (5.2a) 

♦,(?)- 1  +  5?  ?>0  (5.2b) 

Equation  (5.1)  is  expected  to  over  predict  our  dissipation  estimates  as  long  as  there  is  a 
net  flux  of  energy  into  the  ocean.  Evidence  for  this  is  shown  in  Figure  5.6,  which 
illustrates  that  the  measured  dissipation  profiles  are  less  than  their  MO-similarity 
prediction.  These  measurements  are  averaged  over  a  2-hour  period  of  increasing  winds 
when  the  atmosphere  was  providing  energy  to  growing  (i.e.,  young)  seas.  They  were 
taken  during  the  Office  of  Naval  Research  Marine  Boundary  Layers  (MBL)  experiment 
aboard  the  R/P  FLIP  (Edson  et  al.,  1997).  The  data  collected  from  the  FLIP  mast  allows 
us  to  compute  all  the  terms  of  the  total  kinetic  energy  flux  at  three  levels  and  the 
dissipation  rate  of  KE  at  four  levels. 

The  effect  of  the  waves  on  the  surface  energy  budget  can  be  quantified  by  examining 

8 

the  ratio  of  the  predicted  to  measured  dissipation  rate,  — ,  as  a  function  of  sea  state.  The 
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Figure  5.6  Dissipation  measurements  as  a  function  of  height  taken  during  the  ONR’s 
Marine  Boundary  Layer  (MBL)  experiment  aboard  the  R/P  FLIP  compared  to  their 
Monin-Obukhov  similarity  prediction.  The  solid  and  dashed  curves  represent  the 
limits  of  the  standard  error  of  the  measured  and  predicted  values,  respectively. 
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8 

energy  input  into  the  ocean  is  investigated  by  bin  averaging  the  ratio  —  by  the  wave-age 

e, 

C 

parameter,  — .  This  allows  examination  of  the  dissipation  profiles  over  developing, 

u. 

fully-developed,  and  decaying  seas.  The  decay  of  the  wave-induced  effects  as  a  function 
of  sea-state  is  accounted  for  by  multiplying  the  height  of  each  measurement  by  the 
wavenumber  of  the  dominant  wave,  k. 

g 

These  bin-average  estimates  of  —  are  plotted  versus  the  non-dimensional  height 

bp 

C 

kz  in  Figure  5.7.  The  lowest  values  of  —  represents  the  youngest  seas.  The  second  and 

Ut 

third  bins  span  those  wave-ages  that  are  normally  associated  with  fully-developed  seas. 
The  remaining  ranges  correspond  to  decaying  conditions.  These  results  indicate  that 
there  is  an  appreciable  difference  between  the  dissipation  estimates  and  their  MO- 
similarity  predictions  for  all  sea-states  when  kz  <  1 .  The  results  indicate  that  the  deficit  is 
greatest  over  the  youngest  seas  and  decreases  with  increasing  wave  age,  a  conclusion  that 
is  also  supported  by  the  drag  coefficient  computations  which  indicated  increasing 
accuracy  of  the  inertial  dissipation  method  as  the  seas  become  fully-developed.  This 
supports  the  hypothesis  that  the  dissipation  deficit  is  a  function  of  sea  state  and  clearly 
shows  that  the  dissipation  rates  do  not  obey  MO-similarity  in  the  wave  boundary  layer. 
The  results  also  indicate  that  MO-similarity  is  again  valid  above  the  wave  boundary  layer 
at  heights  corresponding  to  kz  >1.5. 
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Figure  5.7  Ratio  of  predicted  to  measured  dissipation  rates  versus  non-dimensional 
height  for  six  wave-age  bins.  Symbols:  o-  -  -  cp/u.  < 24 ;  O —  24 <cp / u .  < 28 

□ - 28 < Cp fu.  < 32 ;  V - 32<Cp/u.  <37 ;  ☆ - 37<Cp/u.  <50;  and  A- - 

50  <  Cp  fu.  <  180 .  The  error  bars  denote  the  standard  error  for  each  measurement. 


Application  of  the  inertial  dissipation  method  to  over  ocean  measurements  therefore 
requires  modification  to  equation  (5.1)  to  account  for  the  wave-induced  effects.  A  simple 
means  to  accomplish  this  involves  fitting  curves  to  the  dimensionless  profiles  shown  in 
Figure  5.7  as 


(5.3) 


where  the  function,  / ,  includes  the  effects  of  the  waves  on  the  production  and  transport 
terms.  Clearly,  the  function  should  decay  to  zero  above  the  wave  boundary  layer.  Its 
behavior  near  the  surface  is  much  less  obvious.  However,  some  insight  into  its  behavior 
is  gained  by  considering  the  surface  energy  flux,  F(o),  which  is  continuous  across  the 
interface.  Therefore,  the  vertical  derivative  of  this  flux,  which  is  closely  related  to  the 
energy  budget  and  e ,  must  approach  zero  very  near  the  surface.  A  simple  function  that 
meets  these  upper  and  lower  constraints  while  closely  fitting  the  measurements  as  shown 
in  Figure  5.8  is  given  by 

/  —  ,kz  =^L(2kz)e~2h  (5.4) 

l  «•  J  Cp 

This  function  is  used  to  correct  the  CMO  dissipation  estimates  in  the  following  analysis 
such  that  the  proposed  dimensionless  dissipation  function  is  now  given  by 

(2kz)e-2kz  (5.5) 

where  <|>E  is  calculated  using  equations  5.2(a)  and  5.2(b). 
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Figure  5.8  Ratio  of  predicted  to  measured  dissipation  rates  versus  non-dimensional 
height  for  six  wave-age  bins.  The  curves  represent  the  function  /  =  ^-(2 Az)r2fc  and 


the  error  bars  denote  the  standard  error.  The  symbols  are  the  same  as  in  Figure  5.7. 


The  modified  inertial  dissipation  estimates  are  compared  to  the  original  dissipation 
estimates  and  the  direct  covariance  measurements  of  the  friction  velocity  in  Figure  5.9. 
Whereas  the  previous  inertial  dissipation  method  systematically  under  predicted  the 
momentum  flux,  the  new  parameterization  improves  the  estimates  substantially.  The 
mean  deviation  of  the  inertial  dissipation  values  calculated  from  the  previous  model  was 
0.049  m/s.  The  estimates  calculated  from  the  proposed  inertial  dissipation  method 
deviate  by  one-third  the  previous  error;  i.e.,  the  mean  deviation  is  only  0.032  m/s.  Thus, 
employing  the  new  parameterization  results  in  over  a  33%  reduction  in  the  error 
associated  with  calculating  friction  velocity  estimates  using  the  inertial  dissipation 
method. 

Importantly,  the  new  parameterization,  <|>En ,  was  developed  using  an  independent 

data  set  measured  over  open-ocean  conditions  of  unlimited  fetch  and  depth.  This 
suggests  the  universality  of  the  modified  function  and  serves  as  satisfactory  validation  of 
the  parameterization.  That  is,  the  parameterization  has  not  been  “tuned”  to  fit  the  CMO 
measurements,  but  rather  was  developed  independently  and  found  to  result  in  impressive 
improvement  of  the  dissipation  method  when  compared  to  the  motion-corrected  direct 
covariance  measurements. 

5.4.  Modified  Bulk  Aerodynamic  Methods 

The  approach  proposed  in  this  section  is  based  upon  the  ideas,  presented  in  section 
4.2.2,  that  Chamock’s  expression  for  the  surface  roughness,  zr ,  accounts  for  the 
roughness  of  the  short  gravity-capillary  waves.  However,  modifications  of  the  drag  due 
to  gravity  waves  and  swell  are  not  included  in  standard  bulk  formula.  To  improve  drag 
coefficient  estimates  and,  therefore,  wind  stress  estimates  in  the  bulk  aerodynamic  model. 
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Inertial  Dissipation  u,  (m/s) 


Friction  Velocity  u*  Comparison 


Figure  5.9  Friction  velocity,  u, ,  comparison  for  original  and  proposed  inertial 

dissipation  methods  and  direct  covariance  method.  Symbols:  o  Proposed  inertial 
dissipation  method;  x  Original  inertial  dissipation  method.  Data  are  bin  averaged 
in  increments  of  0.02  m/s. 
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we  need  to  relate  the  total  drag  (short  plus  long  waves)  to  wave  parameters  in  a 
quantitative  manner. 

The  goal  of  directly  relating  the  sea  surface  drag  coefficient  to  sea  state  has  remained 
elusive  (Smith  et  ah,  1992;  Dobson  et  al.,  1994).  Donelan  (1982)  was  one  of  the  first  to 
present  field  data  which  showed  a  strong  relationship  between  wave  parameters  and  wind 
stress.  His  research  demonstrated  the  breakdown  of  Chamock’s  relation  between  the 
aerodynamic  roughness  length  and  the  friction  velocity  in  all  but  fully-developed  wave 
conditions. 

Since  this  pioneering  work,  other  researchers  have  attempted  to  define  a  quantitative 
relationship  between  the  wind  and  wave  fields.  Smith  (1991)  demonstrated  a  “drag 
coefficient  anomaly,”  a  departure  of  the  drag  coefficient  from  values  expected  through  the 
use  of  Chamock’s  expression,  which  was  significantly  correlated  with  wave  age: 


10  ACOJV  =1.85-2.24 


(ul0N  cose) 


(5.5) 


where  Cp  and  Ul0N  are  expressed  in  ms'1  and  0  is  the  difference  between  the  wind  and 

wave  directions.  Geemaert  et  al.  (1986)  proposed  an  empirical  relation  between  wind 
stress  and  wave  age;  however,  lacking  wave  spectral  data,  they  were  forced  to  estimate 
the  wave  field  from  wind  data  and  may  have  predetermined  their  results  on  account  of 
this  circular  argument.  Following  the  ideas  of  Kitaigorodskii  (1973)  and  Janssen  (1989), 
Nordeng  (1991)  proposed  a  “wave  age  dependent  Chamock  constant”: 


a  =z  =^ 

uc  *■*  2 

u. 


=f\ 


(C. 


yu.  j 


(5.6) 
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However,  this  approach  cannot  convincingly  stand  alone  since  u,  appears  in  both  z,  and 


the  wave  age,  and  self-correlation  may  produce  spurious  results  (Hicks,  1978;  Dobson 
etal,  1994). 

5.4.1.  Scaling  Roughness  Length  with  Wave  Amplitude 

An  alternate  approach,  which  is  explored  here,  is  derived  from  Donelan  et  al.  (1993). 
The  aerodynamic  roughness  length  is  scaled  by  the  rms  amplitude  a  of  the  waves,  where 

z 

a  is  defined  in  terms  of  the  significant  wave  height,  H,  as  4c?  =  H  and  where  — 

a 

expresses  the  ability  of  the  waves  to  serve  as  roughness  elements  (Donelan,  1990).  The 

Z 

benefit  of  regressing  the  ratio  of  —  to  inverse  wave  age  is  that  the  dimensionless  ratios 

CT 

are  formed  from  four  independent  parameters,  and  the  regression  avoids  the  problem  of 
self-correlation  by  estimating  the  sea  surface  roughness  from  wave  age  parameters  (cr 
and  Cp )  that  can  be  measured  in  situ  (Smith  et  al.,  1992). 

Since  the  idea  that  the  roughness  length  should  be  correlated  with  the  height  of  the 
waves  depends  on  the  concept  of  rough  aerodynamic  flow,  attention  is  limited  to  cases  for 
which  Uw  >  7.5  m/s  (Wu,  1980;  Donelan,  1990).  Thus,  the  following  roughness  length 
determination  applies  only  to  fully  rough  flow.  This  formulation  does  not  apply  to  cases 
for  which  the  surface  layer  appears  to  be  aerodynamically  smooth,  Ul0  <  3  m/s,  or  to 

transitional  flow. 

As  depicted  in  Figure  5.10,  inverse  wave-age  estimates  from  the  modified  inertial 
dissipation  method  exhibit  an  increase  of  roughness  scaled  with  rms  wave  height,  as 
expressed  quantitatively  by  the  regression  line 
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(5.7) 


which  displays  good  agreement  to  the  regression  line 


—  =  5.5  xlO”4 
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12.7 
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a 


(5.8) 


fitted  to  data  obtained  in  12  m  of  water  on  Lake  Ontario  (Donelan  et  al.,  1993)  and  agrees 
well  to  the  \uwN  jCp  f  5  ]  power  law  obtained  by  Smith  et  al.  (1992)  from  data  obtained 
during  the  Humidity  Exchange  over  the  Sea  (HEXOS)  experiment  in  1 8  m  of  water  over 
the  North  Sea.  In  contrast,  Dobson  et  al.  (1994)  found  a  weaker  ^Ul0N/Cp  ) 7  ]  power 


law  for  data  obtained  during  the  Grand  Banks  ERS-1  SAR  Wave  Validation  Experiment. 
However,  the  scatter  in  their  data  set  did  not  allow  them  to  resolve  the  slope  (power  law) 
well,  and  they  indicate  that  a  neutral  or  two  way  regression  would  give  a  steeper 
[(c/,0 w/Cp  f  5  ]  power  law,  which  would  be  in  better  agreement  with  the  power  laws 

obtained  using  the  CMO,  Lake  Ontario,  and  North  Sea  data  sets.  The  close  agreement  in 
the  power  laws  obtained  in  four  independent  experiments,  which  ranged  in  environmental 
conditions  from  lake  to  open  ocean,  suggests  the  strength  of  the  relationship  between  the 
surface  roughness  length,  the  rms  height  of  the  dominant  seas  and  the  wave  age,  as 
expressed  in  equation  (5.7). 

In  Figure  5.1 1,  the  CMO  rms  wave  height  measurements  are  empirically  related  to 
the  wave  age  and  wind  speed  by  the  regression  line 

r  ”1-1.8 


_22  T  T 

8  =  2.9xl0~3  — 

(U  cos0)4  Cp 


(5.9) 
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Regression  based  on  Modified  Inertial  Dissipation  Estimates 


Figure  5.10  The  ratio  of  measured  roughness  length  zr  to  rms  wave  height  a  versus 
inverse  wave  age  Ul0  ICp ,  binned  in  increments  of  log(0.035) .  Symbols: .  CMO 

data;  —  regression  on  present  data,  equation  (5.7);  —  Dobson  et  al.  (1994); 

- Donelan  et  al.  (1993),  equation  (5.8); - Smith  et  al.  (1992). 


Regression  based  on  Modified  Inertial  Dissipation  Estimates 


u,n/c 

10  p 


Figure  5.11  The  ratio  of  the  square  of  rms  wave  height  a2  times  the  square  of  gravity 
g 2  versus  inverse  wave  age  UWN  ICp ,  binned  in  increments  of  log10 (0.035) . 

Symbols:  •  •  CMO  data;  - regression  on  present  data,  equation  (5.9); 

- Donelan  et  al.  (1985). 


where  g  is  the  acceleration  of  gravity  and  0  is  the  angle  between  the  winds  and  the 
waves.  Donelan  et  al.  (1985)  found  a  steeper  [(f/10A,  / Cp  )”3  3  ]  relationship  using  their 

Lake  Ontario  field  data.  Examination  of  Figure  5.1 1  reveals  that  the  regression  used  to 
obtain  equation  (5.9)  is  a  relatively  tight  fit.  The  marked  difference  between  the  power 
laws  derived  from  the  Lake  Ontario  and  CMO  data  sets  is  a  reflection  of  the  different 
environmental  regimes,  the  data  from  Donelan  et  al.  (1985)  representing  shallow  lake 
conditions  and  our  data  representing  ocean  conditions  over  the  continental  shelf.  These 
differences  can  include,  for  example,  a  dependence  on  fetch  that  is  felt  at  the  lake  site  but 
not  at  the  CMO  site. 

Assuming  that  the  wind  and  wave  fields  are  aligned  at  the  surface  (cos  9  « l) , 
equation  (5.9)  may  be  simplified 

-  -.1-  -1-0.92 

o  =  5.4xl0“2  —  (5.10) 

.  s  J|_  CP 


Substituting  the  expression  for  the  rms  wave  height,  equation  (5.1 1),  into  equation 
(5.7)  provides  an  expression  for  the  roughness  length 


(5.11) 


which  has  a  steeper  power  law  than  the  expression  derived  by  Donelan  et  al.  (1993) 


(5.12) 


due  to  the  different  empirical  relationship  for  the  dependence  of  the  wave  height  on  wind 
speed.  Equation  (5.11)  implies  that  the  roughness  length  depends  quadratically  on  the 


62 


wind  and  wave-age.  However,  as  mentioned,  this  relationship  for  the  roughness  length  is 
only  applicable  to  the  same  type  of  ocean  conditions  present  at  the  CMO  site  because  the 
parameterizations  do  not  give  the  same  rms  wave  height,  a ,  for  the  same  wind  speed  and 
wave  age.  This  severely  limits  the  widespread  application  of  this  type  of 
parameterization.  Thus,  it  is  preferable  to  return  to  equation  (5.7)  which  appears  to  be  a 
more  robust  relationship  since  the  scaling  gives  similar  results  under  similar  conditions  at 
all  sites. 

Equation  (5.7)  is  now  used  in  the  bulk  aerodynamic  method  as  an  alternate  means  of 
calculating  the  roughness  length.  It  was  hypothesized  that  including  sea  state  dependence 
into  drag  coefficient  calculations  would  improve  the  wind  stress  estimates,  particularly  in 
the  case  of  young  (developing)  seas.  Friction  velocity  estimates  calculated  by  the 
modified  bulk  method  are  depicted  in  Figure  5.12.  Overall,  the  new  parameterization 
slightly  improves  the  u,  estimates  when  compared  to  their  former  values,  although  the 
difference  is  not  statistically  significant.  Of  note,  the  modified  bulk  method  results  in 
higher  friction  velocity  estimates  than  the  original  method  for  all  wind  speeds 
u,  >  0.2  m/s,  and  the  new  parameterization  produces  estimates  that  over  predict  the  actual 
momentum  flux  for  u,  >  0.3  m/s.  Although  statistically  the  modified  bulk  aerodynamic 
method  yields  the  same  error  as  the  original  method,  the  trend  is  significantly  different. 
Whereas  the  original  method  tended  to  over  predict  the  momentum  flux  over  decaying 
seas  and  underestimate  the  flux  over  developing  seas,  the  modified  method  tends  to  over 
predict  the  flux  for  all  wave  ages,  which  suggest  that  there  is  an  inherent  bias  in  the 
correction. 
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Bulk  Aerodynamic  u*  (m/s) 


Friction  Velocity  u*  Comparison 


Figure  5.12  Friction  velocity,  u, ,  comparisons  of  the  bulk  aerodynamic  versus  motion- 
corrected  covariance  method  for  the  original  bulk  method  (+)  and  the  modified 
parameterization  (o).  Values  of  u,  are  bin-averaged  in  increments  of  0.02  m/s. 
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One  limitation  of  this  model  is  the  log-linear  regression  required  by  equation  (5.7). 

A  closer  inspection  of  Figure  5.11  reveals  that  the  points  indicate  that  either  a  curve  or 

two  different  linear  regressions  would  more  adequately  estimate  — .  Additionally,  notice 

CT 


that  the  young  seas. 


U, 


ION 


*1 

vc,  , 


,  only  represent  about  10%  of  the  total  data,  yet  tend  to 


increase  the  slope  (power  law)  significantly.  In  fact,  when  only  the  young  data  is 
considered,  the  power  law  increases  to  [(t/10lV/Cp  )~3  0  ],  whereas  if  only  the  fully- 


developed  and  decaying  seas  are  considered,  the  slope  decreases  to  [(f/10jV/Cp  )"2  0  ]. 

Notice  that  these  results  are  within  the  span  of  ranges  found  by  other  researchers  (Smith 
et  al.,  1992;  Donelan  et  al.,  1993;;  Dobson  et  al.,  1994)  and  could,  in  part,  explain  the 
differences  that  are  found. 

In  Figure  5.13,  the  new  bulk  aerodynamic  estimates  are  compared  to  direct 
covariance  drag  coefficient  measurements.  In  contrast  to  the  bulk  estimates  in  Figure  5.4, 
the  modified  roughness  length  calculation  results  in  sea  state  dependent  drag  coefficients; 
for  a  given  wind  speed,  the  predicted  drag  coefficient  is  inversely  related  to  wave-age. 
However,  the  tendency  of  the  ^l0N/Cp  )'2  7  ]  power  law  to  overestimate  the  friction 

velocity  results  in  higher  than  expected  drag  coefficients  for  all  wave  ages  and  wind 
speeds. 


5.4.2.  Modified  Charnock’s  Relation 


The  overall  goal  is  to  provide  the  best  momentum  flux  estimates  for  the  Coastal 
Mixing  and  Optics  experiment.  While  the  results  obtained  using  the  approach  of  Donelan 
et  al.  (1993)  are  promising,  the  new  estimates  are  not  significantly  improved  over  the 
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tOl-X  D  cOt x 


Drag  Coefficient  vs  Wind  Speed  Comparison 


Cp/u,<=24 


+  + 
o  o 


24<Cp/u,<=28 


32<Cp/u, 


Bulk  Aerodynamic 
Direct  Covariance 


Figure  5.13  A  comparison  of  the  neutral  drag  coefficient  CDN  estimated  using  the 
modified  bulk  aerodynamic  method  and  compared  to  motion-corrected  covariance 

Cm 

measurements  (averaged  in  increments  of  1  m/s)  for  four  different  wave-age  -£■  bins. 
The  error  bars  for  the  covariance  measurements  denote  the  standard  error. 
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previous  values.  Furthermore,  the  parameterization  tends  to  overestimate  the  roughness 
length  for  all  wind  speeds,  a  trend  which  indicates  an  inherent  bias  in  the  model.  Thus, 
still  in  pursuit  of  a  means  to  reconcile  the  current  deficiency  of  the  bulk  aerodynamic 
method  to  estimate  the  wind  stress  in  cases  of  young  (developing)  seas,  but  not  to 
overestimate  the  wind  stress  in  cases  of  fully-developed  or  decaying  seas,  another  method 
of  forcing  a  wave-age  dependence  into  the  roughness  length  via  the  Chamock’s 
“constant”  is  explored.  On  the  basis  of  the  direct  covariance  measurements,  the 
Chamock’s  “constants”  that  provided  the  best  fits  for  young,  fully-developed,  and 
decaying  seas  are  ac  =0.0150,  ac  =0.0126  and  ac.  =  0.0059 ,  respectively.  Figure  5.14 

shows  the  measured  Chamock’s  parameter  as  a  function  of  wave  age,  for  which  the 
function 


a„  =-4.05x10' 


+  2.34x10 


-2 


(5.13) 


provides  the  best  linear  regression.  Revisiting  equation  (4.19),  the  approach  of  Smith 
(1988)  is  modified 

z0=zs+zr  (5.14) 


where  zs  is  determined  as  before  [equation  (4.20)],  but  zr  is  now  expressed  as 

2 


W. 


where 


a„  = -4.05xl0~ 


cA 


\u-  j 


+  2.34x10 


-2 


C„ 

for  —  <  43 

u. 


(5.15) 
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F igure  5.14  Estimated  wave  age  dependent  Chamock’ s  “constant” : 

(C  \ 

"  +  2.34x  1(T2 . 


a,  =  -4.05x10' 


V  u*  J 
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For  seas  with  wave  ages  older  than  —  >  43 ,  the  Chamock’s  parameter  is  assigned  the 

u. 

value  of  ac  =  0.0059 ,  the  value  of  the  upper  limit. 

The  modified  Chamock’s  parameterization  results  in  improved  friction  velocity 
estimates,  as  illustrated  in  Figure  5.15.  Whereas  the  original  bulk  estimates  had  a  mean 
deviation  of  0.046  m/s,  the  new  parameterization  reduces  the  error  by  7%  and  results  in  a 
mean  deviation  of  0.043  m/s.  While  this  may  seem  small,  the  subsequent  improvement  in 
the  drag  coefficient  estimates  (Figure  5.16)  is  quite  noticeable.  Previously,  the  bulk 
aerodynamic  method  underestimated  the  drag  associated  with  young  seas,  but  now  the 
estimates  are  in  good  agreement  with  the  direct  covariance  measurements.  For  the 
remaining  three  wave-age  bins,  the  modified  Chamock’s  expression  again  provides 
improved  estimates  of  the  drag  over  fully-developed  and  decaying  seas.  These  results 
suggest  that  there  is  a  benefit  to  incorporating  wave-age  dependence  into  the  bulk 
aerodynamic  model,  even  if  the  forcing  is  as  simple  as  including  a  wave-age  dependent 
Chamock’s  parameter.  Obviously,  further  research  is  required  to  validate  the  function 
proposed  here,  especially  since  the  wave-age  span  is  limited:  the  youngest  seas  averaged 

Cp  Cp 

~  20 ,  and  the  oldest  seas  averaged  — -  «  43 ,  providing  a  somewhat  limited  range 
w.  u, 

upon  which  the  proposed  function  is  based.  Thus,  while  the  relationship  proposed  here 
improves  the  bulk  aerodynamic  momentum  flux  estimates  for  the  CMO  experiment, 
future  studies  are  needed  to  validate  or  improve  the  parameterization.  Most  importantly, 
the  results  indicate  that  bulk  aerodynamic  momentum  flux  estimates  can  be  improved  by 
incorporating  sea  state  dependence  via  a  roughness  length  determination. 
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Bulk  Aerodynamic  u,  (m/s) 


Friction  Velocity  u*  Comparison 


Figure  5.15  Friction  velocity,  ut ,  comparisons  of  the  bulk  aerodynamic  versus  motion- 
corrected  covariance  method  for  the  original  bulk  method  (+)  and  the  modified- 
Chamock  parameterization  (o).  Values  of  ut  are  bin-averaged  in  increments  of 
0.02  m/s. 
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Drag  Coefficient  vs  Wind  Speed  Comparison 


Cp/u.<=24  24<Cp/u„<=28 


+  +  Bulk  Aerodynamic 

o  o  Direct  Covariance 


Figure  5.16  A  comparison  of  the  neutral  drag  coefficient  CDN  estimated  using  the 
modified-Chamock  bulk  aerodynamic  method  and  compared  to  motion-corrected 
covariance  measurements  (averaged  in  increments  of  1  m/s)  for  four  different  wave- 

age  bins.  The  error  bars  for  the  covariance  measurements  denote  the  standard 
error. 


71 


6.  Conclusions 


The  work  presented  in  this  thesis  is  an  investigation  of  the  influence  of  surface  waves 
on  momentum  exchange.  Accurate  measurements  of  the  exchange  of  momentum  and 
energy  between  the  atmosphere  and  ocean  are  critical  to  numerical  prediction  systems. 
The  growing  sophistication  of  current  numerical  models,  such  as  the  United  States 
Navy’s  Coupled  Ocean  Atmosphere  Prediction  System,  requires  an  accurate 
representation  of  the  dynamics  and  magnitude  of  energy  transfer  through  the  air-sea 
interface.  This  thesis  explored  the  theoretical  framework  which  forms  the  basis  for  the 
direct  covariance  (eddy  correlation),  direct  /  inertial  dissipation  and  bulk  aerodynamic 
methods  for  measuring  and  calculating  momentum  flux.  Unlike  the  other  methods  which 
rely  on  similarity  theory,  the  direct  covariance  method  provides  a  direct  computation  of 
the  wind  stress.  However,  this  method  is  highly  susceptible  to  flow  distortion  and 
platform  motion  and  is,  therefore,  difficult  to  implement  in  over-ocean  studies.  It  is 
therefore  preferable  to  develop  indirect  methods  such  as  the  bulk  aerodynamic  and 
inertial  dissipation  methods  which  relate  easily  measured  quantities  such  as  wind  speed 
and  temperature  to  the  momentum  flux.  These  methods  allow  us  to  estimate  turbulent 
quantities  when  direct  measurements  are  unavailable  and  to  estimate  the  lower  boundary 
flux  conditions  to  input  into  numerical  models. 

The  indirect  methods  for  estimating  momentum  exchange  rely  on  Monin-Obukhov 
similarity  theory,  which  was  developed  almost  half  a  century  ago  and  has  been  validated 
over  land  by  a  number  of  field  experiments.  Its  applicability  over  the  ocean  has  only 
recently  been  examined.  By  quantitatively  comparing  the  bulk  and  inertial  dissipation 
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estimates  to  direct  measurements,  it  was  determined  that  the  indirect  methods 


systematically  underestimate  the  momentum  flux  into  developing  seas.  It  was 
hypothesized  that  incorporating  sea  state  dependence  into  the  parameterizations  would 
improve  the  wind  stress  estimations  computed  from  both  methods. 

For  the  inertial  dissipation  method,  a  new  parameterization  for  the  dimensionless 
dissipation  rate  that  includes  the  effects  of  the  surface  waves  and  accounts  for  the  non¬ 
zero  flux  of  energy  from  the  atmosphere  to  the  ocean  is  proposed: 


This  parameterization  was  developed  from  an  independent  data  set  and  has  resulted  in 
over  a  33%  reduction  of  the  error  associated  with  the  inertial  dissipation  friction  velocity 
estimates  when  compared  to  direct  covariance  measurements. 

Two  separate  methods  of  incorporating  sea  state  dependence  into  the  bulk 
aerodynamic  method  were  explored.  The  first  parameterization  is  based  on  the  ideas  of 
Donelan  et  al.  (1993)  and  involved  the  derivation  of  a  sea  surface  roughness  length 
normalized  by  rms  wave  height  that  is  a  log-linear  function  of  wave  age: 


(6.2) 


This  parameterization  agrees  well  with  the  relationships  presented  by  Smith  et  al.  (1992) 
and  Donelan  et  al.  (1993). 

The  parameterization  was  used  in  the  bulk  aerodynamic  method  and  slightly 
improves  the  friction  velocity  estimates,  although  the  difference  is  not  statistically 
significant.  The  results  indicate,  however,  an  inherent  bias  and  systematic  trend  in 
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overestimating  w, .  One  limitation  is  the  lack  of  data  with  young  wave  ages  so  that  the 
power  law  cannot  be  resolved  well  in  the  lower  limit.  Another  possibility  is  that  the  log- 
linear  relationship  is  too  simplistic;  the  data  suggest  that  a  curve,  or  at  least  two  lines, 
would  provide  a  more  optimal  fit. 

The  second  approach  to  modifying  the  bulk  aerodynamic  method  involved 
incorporating  a  sea  state  dependent  Chamock’s  parameter  into  the  roughness  length 
calculation: 


z 


r 


(6.3a) 


where 


a„  =  -4.05x10' 


V^lOAf 


+  2.34x10 


-2 


(6.3b) 


C 

for  —  <  43 . 

m. 


When  compared  to  the  direct  covariance  measurements,  this  new  parameterization 
resulted  in  a  7%  reduction  in  the  error  associated  with  bulk  aerodynamic  friction  velocity 
estimates  and  considerably  improved  drag  coefficient  estimates.  The  results  are  tentative 
due  to  the  limited  wave-age  span  upon  which  the  function  is  derived,  and  more  research 
is  needed  to  validate  or  improve  the  shape  of  the  function,  especially  in  both  the  upper 
and  lower  limits. 

Finally,  this  research  strongly  supports  the  use  of  either  the  bulk  aerodynamic  or 
inertial  dissipation  methods  to  estimate  momentum  exchange  over  the  ocean.  Since 
widespread  direct  measurements  are  not  practical,  indirect  methods  such  as  the  bulk  and 
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inertial  dissipation  methods  are  invaluable.  While  future  research  will  undoubtedly 
improve  the  accuracy  of  both  methods,  especially  as  we  gain  increased  understanding  of 
wind- wave  coupling,  the  results  obtained  here  indicate  that  the  current  bulk  aerodynamic 
and  inertial  dissipation  algorithms  which  already  have  widespread  acceptance  are,  in  fact, 
doing  an  adequate  job  of  estimating  the  momentum  flux.  This  investigation  provides 
strong  quantitative  evidence  that  incorporating  sea  state  dependence  into  both  methods 
takes  us  one  step  closer  and  improves  upon  current  momentum  flux  estimates. 
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